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Abstract 


A  model  is  developed  which  treats  the  growth  o-f  a 
'fireball  from  an  atmospheric  nuclear  detonation  up  to  the 
point  of  hydrodynamic  separation.  The  radius  of  the  fireball 
and  the  partition  of  energy  among  streaming  radiation, 
material  temperature,  and  thermal  radiation  are  predicted  as 
a  function  of  time.  A  simplified  model  is  developed  which 
applies  when  equilibrium  between  material  temperature  and 
thermal  radiation  temperature  is  reached. 

Numerical  results  are  calculated  for  various  yields  and 
yield  pulse  shapes.  These  results  are  compared  to  those 
predicted  by  the  model  previously  developed  by  6.  C. 
Pomraning. 

The  model  developed  here  predicts  reasonable  numerical 
results  and  physically  reasonable  time  dependences.  It  does 
not  entail  integrations  over  retarded  time  or  the  use  of 
moments  of  radiation  transport  equations.  It  engenders  an 
intuitive  understanding  of  the  fireball  growth  process  and 
was  developed  for  pedagogical  use. 


1  ■  Introductiixi 


A.  Background 

1.  Fireball  Radi ativa  Growth  Phaaa.  This  thesis  deals 
with  the  radiative  growth  phase  of  the  developsent  of  the  m- 
ray  fireball  from  an  atmospheric  nuclear  explosion.  During 
this  phase,  x-rays  from  the  device  heat  the  surrounding  air 
to  a  high  level  of  ionization.  A  thermal  radiation  field 
builds  up  as  a  result  of  re-emission  of  photons  by  the  air. 
The  fireball  grows  as  photons  from  the  device  stream  outward 
to  be  absorbed  in  the  ambient  air  surrounding  the  fireball. 
The  other  mechanism  of  transport  of  energy  outward  is  the 
diffusion  of  the  thermal  radiation  in  the  fireball.  During 
this  phase,  energy  transport  to  the  surface  of  the  fireball 
is  accomplished  only  by  radiative  transport  because  the 
growth  is  at  a  rate  faster  than  a  material  shock  wave  can 
propagate  through  the  fireball,  even  under  the  extreme 
pressures  and  temperatures.  This  radiative  phase  ends  at 
hydrodynamic  separation,  when  the  radiative  transfer  has 
decreased  enough  that  the  shock  front  can  travel  as  fast  as 
the  fireball  radius. 

Modeling  the  Radiative  Growth.  The  standard 
reference  in  the  field,  "The  Effects  of  Nuclear  Weapons" 
CRef  13  provides  qualitative  information  about  the  radiative 
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growth  o-f  •  nuclear  -fireball,  but  no  natheoatical  model  of 
the  processes.  G.  C.  Pomraning  has  published  such  a 
mathematical  model  CRe-f  23,  referred  to  here  as  the  "early 
time  model",  following  the  title  of  his  paper.  Reference  3 
develops  the  model  in  more  detail,  including  presentation  of 
a  finite  differencing  scheme  for  implementing  the  model.  His 
model  allows  the  computation  of  the  material  temperature  and 
thermal  radiation  temperature  but  is  limited  in  its 
effectiveness  for  pedagogical  use  by  its  complexity,  as  well 
as  for  other  reasons  examined  in  detail  later. 

3.  Pedagogical  Model .  There  exists  the  need  for  a 
mathematically  simple  model  udiich  is  based  on  the  physical 
processes  of  the  radiative  fireball  growth.  This  model  is 
needed  for  pedagogical  use  in  programs  such  as  the  Nuclear 
Weapons  Effects  curriculum  at  the  Air  Force  Institute  of 
Technology,  where  this  research  was  done.  Such  a  model 
should  predict  numerically  and  physically  reasonable  results 
and  should  engender  an  intuitive  physical  understanding  of 
the  nuclear  fireball  radiative  growth. 

Problem  and  Scope 

The  objective  of  this  research  has  been  to  develop  a 
physically  meaningful  model  of  the  growth  of  a  fireball  from 
an  atmospheric  nuclear  detonation  up  to  the  point  of  hydro- 
dynamic  separation  which: 


r 

—  is  appropriats  -for  instruction  at  ths  graduats 
Isvsl  Nithout  rsquiring  a  prsvious  knoMlsdgs  of  transport 
thsory  or  advancsd  sathMatical  or  numsrical  Mithods 

—  is  bassd  on  simpls  sodols  of  ths  physical  procsssss 

—  prsdicts  tsflipsraturss  and  tins  dspsndsncss  of  ths 
right  order  of  magnitude 

—  predicts  curves  of  radius  and  temperatures  Mith 
physically  reasonable  shapes. 

The  second  objective  of  this  research  has  been  to 
verify  these  characteristics  by: 

*—  implementing  the  model  in  a  computer  program 

—  graphing  the  results  for  various  combinations  of 
yield,  pulse  duration,  and  pulse  shape 

comparing  to  results  obtained  using  the  early  time 
model  as  a  benchmark. 

The  scope  of  the  project  has  been  limited  to  the 
consideration  of  the  radiative  transport  phase  only.  The 
examples  computed  for  evaluation  and  comparison  include  t»«o 
nuclear  pulse  shapes,  a  square  eave  and  an  isosceles 
triangle.  The  eKamples  all  assume  sea  level  ambient  air 
density. 

C.  Assumptions 

All  the  models  presented  here  assume  that  the  nuclear 
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dstonation  is  a  -free  air  burst  in  uniform  air,  that  thm 
nuclmar  dsvics  is  a  point  sourcm  of  mnmrgy  in  thm 
atmosphmrm»  and  that  thm  fireball  is  homogeneous. 

Approach  and  Presentation 

The  sequence  of  presentation  of  this  thesis  follows  the 
approach  used  in  the  research.  Section  XI  summarizes  the 
early  time  model  and  presents  the  results  of  its  application 
to  an  example  problem.  Difficulties  with  the  model  are 
considered.  Section  III  considers  the  approach  to  be  used  in 
developing  a  new  model.  Similarities  with  the  early  time 
model,  and  assumptions  in  common  are  presented,  followed  by 
tferences  and  additional  assumptions.  Section  IV  develops 
the  homogeneous  three-temperature  model  using  energy  balance 
equations.  This  model  applies  when  the  temperatures  are  not 
yet  in  equilibrium.  A  finite  differencing  algorithm  is 
presented.  Section  V  develops  the  homogeneous  one- 
temperature  model  which  applies  after  equilibrium  has  been 
established  between  the  material  temperature  and  the  thermal 
radiation  temperature.  Section  VI  presents  the  results  of 
numerical  calculations  using  the  combined  homogeneous  model, 
and  comparing  the  results  to  those  obtained  using  the  early 
time  model.  Section  VII  evaluates  the  pedagogical  and 
numerical  advantages  and  limitations  of  the  homogeneous 
model,  and  conclusions  and  recommendations  are  presented  in 


section  VI I I 


II-  Background;  Po«TanlnQ*»  Early  Tima  Modal 


A.Charactariaticm  o-f  tha  Radiativa  Browth  Phaaa 

Tha  radiativa  groath  phaaa  of  firaball  foraation  could 
ba  dividad  into  subphaaaa  baaad  on  aavaral  posaibla 
distinctions.  Among  thasa  ara  (1)  burnout  vs.  not  fully 
ionizadl  (2)  anargy  transport  aachanisml  (3)  aquilibriua  vs. 
nonaqui librium  of  matarial  tamparatura  and  tharmal  radiation 
taaparatura.  Tha  assumptions  chosan  in  ragard  to  thasa 
char actari sties  largaly  datarmina  tha  natura  of  tha  modal 
darivad.  Tharafora,  thasa  distinctions  will  ba  dascribad 
bafora  tha  assumptions  for  tha  aarly  tima  OMsdal  and  for  tha 
homoganaous  modal s  ara  prasantad. 

I-  Burnout.  At  tha  baginning  of  firaball  formation* 
tha  voluma  of  tha  firaball  is  small.  Tharafora*  tha 
radiation  is  so  intansa  as  to  fully  ioniza*  i.a.  burnout* 
tha  air  in  tha  firaball.  Any  tima  an  alactron  raattachas  to 
a  nuclaus  (amitting  a  photon  of  tharmal  radiation)*  it 
immadiataly  absorbs  anothar  photon  of  straaming  radiation 
and  tha  atom  is  again  fully  ionizad.  Thus*  tha  matarial 
tamparatura*  T^*  is  maintainad  at  tha  burnout  tamparatura* 
T|,»  approx imataly  400  aV  CRaf  23.  As  tha  firaball  bacomas 
largar*  and  (avantually)  as  tha  production  of  straaming 
radiation  by  tha  nuclaar  axplosiva  dacraasas*  tha  radiation 
bacomas  lass  intansa.  At  soma  point*  atoms  which  form  ara  no 
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lang«r  iMaadiataly  r»-ic3niz«d,  and  tha  firaball  cools.  At 
this  point,  tha  burnout  subphasa  is  coaplatad. 

Enargy  Transport  Machanisa.  Thara  ara  two  anargy 
transport  aachanisas  availabla  during  tha  radiativa  growth 
phasa.  Thasa  ara:  straaaing  of  radiation  froa  tha  nuclaar 
axplosiva  diractly  outward  through  tha  firaball,  and 
diffusion  of  tha  tharaal  radiation  producad  in  tha  voluaa  of 
tha  firaball.  At  tha  baginning  of  firaball  foraation,  tha 
straaaing  radiation  is  intansa  and  tha  tharaal  radiation  is 
nagliglbla,  so  tha  aiachanisa  is  straaaing  only.  Evantually, 
tha  yiald  pulsa  is  ovar  and  tha  straaaing  radiation  has  baan 
absorbad,  but  tha  tharaal  radiation  has  built  up 
significantly,  so  tha  aachanlsa  is  diffusion  only.  Tha 
transition  froa  straaaing  only  to  diffusion  only  is  a  saooth 
ona,  so  for  soaa  intarmadiata  intar val  tha  transport 
aachanlsa  is  both  straaaing  and  diffusion. 

Eciullibrlua.  At  tha  baginning  of  firaball 
foraation,  tha  aatarial  is  haatad  (assantially  instantly)  to 
burnout  taaparatura  but  tha  tharaal  radiation  fiald  takas 
tiaa  to  build  up.  Ccxiditions  start  far  froa  aquilibriua: 
T^  ■  Tjj  but  T^y  y  ■  O  ,  whara  T^^  is  tha  affactiva  (black 
body  aquivalant)  taaparatura  of  tha  tharaal  radiation.  Nhila 
tha  straaaing  radiation  is  intansa,  and  tha  firaball  is 
rapidly  growing  (diluting  tha  tharaal  radiation  as  it  builds 
up),  this  nonaquilibriua  condition  will  ba  aaintalnad. 
Evantually,  howavar,  aquilibriua  will  ba  closaly  approachad. 


B.  Major  Aasumptionm  g£  tha  Early  Tima  Modal 

Tha  aarly  tima  modal  dividaa  tha  radiativa  growth  phaaa 
into  two  diatinct  subphasea,  burnout  and  diffuaion.  During 
tha  burnout  aubphaaa,  the  machaniam  of  energy  tranaport  ia 
aaaumed  to  ba  atreaming  only.  During  the  poat -burnout 
pariodf  tha  mechanism  is  aaaumed  to  ba  diffuaion  only. 

The  early  time  modal  aaaumaa  that  atreaming  radiation  la 
abaorbad  at  whatever  rate  ia  required  to  maintain  at 
Thia  ia  equivalent  to  assuming  that  the  air  is  opaque  to 
streaming  radiation  <  fC  ^  ~00)  if  ^  and  transparent  to 
streaming  radiation  < /C .  0)  if  >  T.  .  Thus,  the  model 

assumes  that  the  breakpoint  between  burnout  and  cooling  and 
the  breakpoint  between  streaming  only  and  diffuaion  only 
occur  simultaneously.  The  model  further  assumes  that  there 
is  no  overlap  period  of  transport  by  simultaneous  atreaming 
and  diffusion. 

The  model  assumes  nonequilibrium  of  temperatures 
throughout  both  subphases.  The  model  predictSp  numericallyp 
the  close  approach  to  ef^uilibriump  but  does  not  require 
equilibrium  as  a  simplifying  assumption. 

The  early  time  model  assumes  that  tha  burnout  phase 
ends  when  one  of  the  following  conditions  is  met  (whichever 
occurs  first) : 

—  The  streaming  radiation  intensity  at  the 
surface  of  the  fireball  falls  to  zerop  or 


All  of  th»  yiold  dolivorod  so  far  is  accountsd 
for  as  Mtsrial  and  tharmal  radiation  onsrgy  of  ths  firsball. 

Tho  sarly  tins  oodsl  usss  vory  different  tschniciusB 
for  aiodsling  ths  two  phases,  burnout  and  diffusion.  Those 
two  approaches  are  discussed  next. 

Burnout  Phase  Approach 

Consider  the  fireball  at  a  tiine  t,  at  which  it  has  a 

radius  R,  with  streaming  radiation  energy  density  E^(r,t). 

Then  energy  which  streams  through  the  boundary  at  R  in  a 

2 

time  increment  dt  is  R  c  E^(R,t)  dt.  This  energy  burns 

out  an  additional  shell  of  thickness  dR  containing  material 

2 

and  streaming  energy  of  4  7f  R  dR  (E^(R,t)  •*>  E^>,  where  E|^ 
is  the  energy  density  of  the  air  at  burnout.  An  energy 
balance  of  these  two  expressions  can  be  solved  for  dR/dt,  if 
E^(R,t)  is  known. 

The  source  of  E^<R,t>  is  yield  which  was  delivered  at  a 
retarded  time  t*  s  t  -  R/c  and  which  is  subject  to 
absorption  as  it  travels  from  r  ■  O  (at  t'>  to  r  ■■  R  (at 
t).  Therefore,  an  integral  over  this  path  is  required  to 
find  the  absorption  loss.  This  path  integral  can  be  treated 
equivalently  as  an  integral  over  retarded  time  from  t*  to  t. 
As  the  streaming  radiation  passes  a  radius  r”  at  time  t" 
(dummy  integration  variables),  the  absorption  is  found  by 
the  requirement  that  it  be  exactly  enough  to  maintain  the 
material  temperature  at  burnout,  i.s>  T_  ■  T.  so  the  energy 

fli  D 


B 


absorbed  from  streaming  radiation  equals  that  added  to 
thermal  radiation.  The  net  production  rate  of  thermal 
radiation  depends  on  the  emission  rate,  which  is  constant 


.t 


since  T  »T.  ,  and  the  absorption  rate,  which  is  proportional 
to  the  thermal  radiation  energy  density,  E^.  8o,  knowledge 
of  E^(t)  is  required  for  the  integration  which  gives 
E^(R,t)y  which  is  required  to  obtain  dR/dt. 

Since  E^  is  produced  at  a  constant  rate  by  emission 
from  the  air  at  it  would  build  up  exponentially  to 
equilibrium,  except  that  it  is  also  subject  to  dilution  as  R 
increases.  So  knowledge  of  E|.<t)  requires  knowledge  of  R(t>. 

The  result  of  the  approach  is  a  model  consisting  of  a 
coupled  system  of  integral  and  differential  equations  which 
must  be  solved  numerically.  Further  details  of  these 
equations  will  not  be  provided  here,  since  this  summary  is 
intended  only  to  provide  a  conceptual  overview  of 
Pomraning’s  work.  Details  can  be  found  in  either  reference  2 
or  reference  3. 

D.  Diffusion  Phase  Approach 

For  the  diffusion  phase,  Pomraning  first  assumes  that 
the  diffusion  equation  is  valid,  and  that  the  material  and 
thermal  energy  densities  are  not  necessarily  spatially 
uniform.  He  further  assumes  that  E^  ■  £  where  £ 
varies  with  time,  but  is  spatially  uniform  at  each  time. 
This  is  combined  with  an  assumed  cubic  temperature 


9 


dapendencB  of  the  mean  free  path  of  thermal  photons  to 
simplify  the  diffusion  equation.  The  first  and  second 
moments  of  the  diffusion  equation  are  then  taken.  These 
moment  equations  are  then  simplified  by  applying  the 
assumption  that  all  the  energy  densities  E.  ,  E_  and  E_ 

W  »  PI 

are  spatially  uniform  inside  and  zero  outside  the  fireball. 

These  assumptions  result  in  a  coupled  system  of 
differential  equations  which  must  be  solv^  nume»rically. 
Again*  the  details  are  left  to  reference  2.  The  results  of 
calculations  using  the  early  time  model  for  an  example  case 
are  presented  next*  followed  by  an  evaluation  of  some 
difficulties  observed. 

Example  Case  Using  the  Early  Time  Model 

As  an  example  of  the  application  of  the  early  time 
model*  a  20  kiloton  burst*  with  a  constant  yield  delivery 
rate  of  3.8873  kilotons/shake  and'  with  S. 14493  shakes 
duration*  was  analyzed.  This  pulse  duration  was  chosen  to 
provide  the  same  yield  delivery  rate  as  the  peak  value  for 
the  double-exponential  pulse  considered  in  reference  2. 

Figures  1  through  3  show  the  material  and  radiation 
temperatures*  the  radius*  and  the  growth  rate*  dR/dt*  over 
the  first  20  shakes*  starting  at  the  beginning  of  the  burst. 
Observations  about  these  curves  are  presented  in  section  F* 
below.  Examples  of  other  yields*  pulse  shapes*  and  pulse 
durations  are  provided  in  Appendix  A. 
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Fig.  1  Material  and  E-f-fectiva  Radiation  Temperature 


Fig.  2  RadiuBf  Early  Time  Model,  20KT 
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Fig.  3  Growth  Rate,  Early  Time  Model,  20KT 

F.  Predictive  Pr9bl.effls  o£  the  Early  Time  Model 

The  exaoiple  curves  presented  above  illustrate  sone 
problems  with  the  results  predicted  by  the  early  time  model. 

E'f'fective  Temperature.  should  be  a  smoothly 

varying  function  of  time,  yet  figure  1  shows  a  sudden  break 
in  slope  at  the  transition  between  burnout  and  diffusion 
phases. 

2.  Decreasing  Fireball  Expansion  Rate.  dR/dt  should  be 
a  monotone  decreasing  function  of  time  unless  the  yield  rate 
is  rapidly  Increasing,  which  is  never  the  case  for  the 
example  square  wave  pulse.  In  particular,  when  the  streaming 
radiation  is  fully  absorbed,  dR/dt  should  decrease  because 
of  the  less  efficient  transport  of  energy  by  diffusion 
alone,  rather  than  diffusion  and  streaming  together. 
However,  figure  3  shows  a  step  increase  in  dR/dt  at  the  end 
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of  the  burnout  phase. 

3.  Constraint  by  the  Speed  of  Light.  The  enpansion 
rate  of  the  fireball  should  never  exceed  c,  the  speed  of 
light  (3  meters/shake).  However,  the  early  time  eodel 
predicts  an  interval  at  the  beginning  of  the  diffusion  phase 
in  which  dR/dt  is  greater  than  the  speed  of  light. 

4.  versus  T  The  material  temperature,  T_, 

“W  ““  Ctt  fli 

should  always  exceed  the  thermal  radiation  effective 
temperature,  since  emission  by  the  material  is  the 

only  source  of  thermal  radiation  and  since  the  radiation  is 
continually  leaking  out  of  the  fireball  interior  whereas  the 
material  energy  cannot  leak.  However,  figure  1  shows  a  rapid 
drop  in  T^  at  the  beginning  of  the  diffusion  phase,  until  T_ 
is  less  than  T^^^.  T^^^  then  approaches  equilibrium  with  T^ 

from  above,  rather  than  from  below. 

6.  Causes  of  the  Predictive  Problems 

1.  Effective  Temperature.  The  sharp  peak  in  T^^ ^  at 
the  end  of  the  burnout  phase  is  caused  by  not  accounting  for 
diffusion  transport  as  becomes  significant  and  E^  becomes 
small.  This  has  two  related  effects:  E^  is  not  lost  (in  the 
model)  by  conversion  to  E_  as  it  diffuses  into  ambient  air 

HI 

and  is  absorbed,  and  dR/dt  is  underestimated,  so  the 
tendency  to  decrease  E^  by  dilution  is  underestimated.  Both 
effects  cause  E^  to  rise  too  steeply  and  too  high  during  the 
latter  part  of  the  burnout  phase. 
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2*  D»cr»aainq  Fireball  Expansion  Rate.  The  step 
increase  in  dR/dt  at  the  end  O'f  the  burnout  phase  arises 
from  this  sane  cause.  dR/dt  is  underestimated  at  the  latter 
part  of  the  burnout  phase  because  the  diffusion  transport  is 
ignored.  The  effect  is  exacerbated  by  the  overesti nation  of 
T^j. y  and  E^t  since  this  leads  to  an  overestimate  of  dR/dt  in 
the  beginning  of  the  diffusion  phase. 

3.  Constraint  by  the  Sgeed  of  Light.  The  early  time 
aiodel  predicts  a  rate  of  expansion  of  the  fireball  ehich 
exceeds  the  speed  of  light  during  the  diffusion  phase 
because  dR/dt  is  not  being  modeled  directly.  The  moments 
equation  method  does  not  reflect  the  physical  constraint  of 
dR/dt  <  c. 

4.  T_  versus  T  ...  The  early  time  model  predicts 
T^  <  T^j. y  during  the  diffusion  phase  because  vLe  moments 
equations  do  not  directly  model  the  process  of  conversion  of 
E.  to  E_  as  the  fireball  expands. 

H.  Fedagogi cal  Difficulties 

1.  Integrals  over  Retarded  Time.  The  retarded  time  and 
path  integrals  can  be  confusing.  The  physics  is  obscured  by 
the  mathematics.  For  the  beginning  graduate  student  with  a 
general  engineering  background*  without  specific  preparation 
in  transport  theory  or  in  the  use  of  retarded  potentials  in 
electromagnetics,  this  can  be  a  problem. 


2.  Moments  of  •  Di-f-fuaion  Equation.  Although  tha 
diffusion  aquation  is  familiar  to  tha  studant,  in  taros  of 
nautron  diffusion  and  Pick's  Law*  tha  usa  of  momants 
aquations  may  be  accepted  (intellectually)  as  a  legitimate 
mathematical  manipulation  but  without  leading  to  an 
intuitive  physical  understanding.  The  moments  equations  are 
particularly  hard  to  accept  when  applied  to  a  homogeneous 
fireball  with  a  step  decrease  to  zero  energy  densities  at 
its  surface*  since  the  diffusion  equation  is  not  valid  near 
the  discontinuity  and  predicts  no  net  energy  transfer 
elsewhere.  A  more  involved  derivation  of  the  same  final 
equations  directly  from  the  full  Boltzmann  equation  might  be 
possible*  but  this  would  present  its  own  problems*  since  the 
student  would  not  have  the  necessary  background  in  transport 
theory*  in  most  cases. 

3.  Computational  Difficulty.  The  early  time  model  does 
not  lend  itself  well  to  use  as  a  calculational  exercise  or 
project.  A  straightforward  finite  differencing  scheme  can  be 
presented*  but  the  integrals  over  retarded  time  require 
maintaining  several  vectors  of  past  variable  values 
amounting  to  an  array  of  several  thousand  floating  point 
numbers.  This  precludes  the  use  of  personal  microcomputers 
and  results  in  more  loading  of  main  frame  facilities. 
Calculation  of  the  numerical  integrations  at  each  time  step 
also  results  in  relatively  long  (eMpensive)  execution  times. 


IS 


Th«  coNiputations  -for  th«  •xamples  in  this  thesis  Mere 
actually  per-formed  on  a  personal  microcomputer,  but  this 
required  the  derivation  o-f  a  system  of  recursive  equations 
which  could  be  iterated  at  each  time  step  to  obtain  the 
needed  retarded  time  integrals,  together  with  the  use  o^ 
programming  tricks  to  increase  speed  at  the  expense  o-f 
intelligibility  o-f  the  program  source  code.  The  recursion 
system  is  presented  in  Appendix  D.  This  would  not  make  the 
model  use-ful  as  a  programming  project,  however,  because  as 
much  time  could  be  spent  presenting  the  programming 
techniques  as  would  be  spent  presenting  the  physics. 

4.  Predictive  Problems.  The  prediction,  numerically, 
of  physically  unreasonable  behaviors,  as  discussed  above, 
presents  problems  in  the  pedagogical  use  of  the  model.  The 
student  would  (one  hopes)  be  unsatisfied  with  an 
understanding  which  is  founded  on  a  model  which  predicts 
unreasonable  behaviors  through  limitations  which  are  not 
understood.  The  difficulty  is  that  to  see  the  causes  of  the 
predictive  problems  of  the  early  time  model  requires  at 
least  as  good  an  understanding  of  the  physics  as  the  model 
is  intended  to  create. 

I .  Conclusions 

A  successful  pedagogical  model  should  avoid  the  use  of 
retarded  time  and  of  moments  equations.  A  finite 
differencing  scheme  should  follow  from  the  model  in  a  direct 
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way.  Tha  di-F-farencing  schema  should  ba  programmabla  as  an 
aMercisa  using  limited  computer  resources,  such  as  personal 
microcomputers.  The  numerical  predictions  of  such  a  program 
should  be  physically  reasonable.  The  modal  should  ba  based 
on  simple  mathematical  models  of  the  relevant  physical 
processes.  The  next  section  will  present  an  approach, 
together  with  its  major  assun^tions,  which  is  then  used  to 
develop  a  model  to  meet  these  goals.  These  assumptions  are 
selected  so  as  to  avoid  those  major  assumptions  which  led  to 
difficulties  with  the  early  time  model. 


III.  Approach  and  AssuAgtigns 


Approach 

Th0  approach  used  in  choosing  assumptions  for  the  model 
to  be  developed  in  subsequent  sections  is  to  retain  as  much 
as  possible  of  the  early  time  model.  Those  assumptions  Mhich 
lead  to  physically  unreasonable  numerical  results  or  to 
pedagogical  difficulties  were  identified  in  the  previous 
section.  Those  assumptions  are  replaced  by  alternative 
choices  as  described  below. 

Assumptions  Retained  from  the  Early  Jime  Model 

Infinite!  Homogeneous  Ambj^ent  Air.  The  energy 
release  from  the  nuclear  explosive  device  is  assumed  to 
occur  in  an  infinite  volume  of  air  of  uniform  thermodynamic 
properties.  This  is  a  reasonable  assumption  because  the 
diameter  of  the  fireball  at  hydrodynamic  separation  is  on 
the  order  of  tens  of  meters  (about  25  meters  for  20  kilo- 
tons)  CRef  1.651  while  the  atmospheric  density  varies 
approximately  exponentially  with  a  scale  height  of  roughly 
7  km  (at  sea  level)  CRef  41.  This  also  assumes  that  the 
burst  is  higher  off  the  ground  than  the  radius  of  the 
fireball  at  hydrodynamic  separation.  The  presence  of  the 
ground  could  be  incorporated  into  the  model,  as  will  be 
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•xplalnsd  in  th»  d*velopn«nt  of  th«  modol,  but  th«  ground  is 
not  sxplicitly  considsrsd  hsrs  bscsuss  it  is  doss  not  add 
significant  physics  to  the  nodel. 

2-  Point  Source  Nuclear  Device.  The  explosive  is 
assumed  to  be  a  point  source  of  energy  release  in  the 
atmosphere.  This  assumption  can  be  viewed  in  either  of  two 
ways.  First,  the  assumption  is  equivalent  to  assuming  that 
the  structure  which  surrounds  the  actual  nuclear  fuel 
behaves  the  same  as  the  air.  Alternatively,  any  effects  of 
this  structure  can  be  considered  to  have  been  incorporated 
in  the  choice  of  the  total  (x-ray)  yield  and  the  normalized 
yield  release  rate  function,  f(t>. 

3.  Homogeneous  Fireball.  The  conditions  inside  the 
fireball,  in  particular  the  energy  densities  and  E,.,  are 
assumed  to  be  spatially  uniform  throughout  the  radiative 
phase  of  the  fireball  growth.  At  burnout,  the  mean  free  path 
of  photons  is  about  50  meters,  which  is  large  compared  to 
the  fireball  radius,  so  the  assumption  is  Initially  an 
excellent  one.  By  the  time  of  hydrodynamic  separation,  at 
about  25  eV  CRef  1.653,  the  mean  free  path  is  about  1/3 
meter,  so  the  homogeneous  assumption  constitutes  a  limita¬ 
tion  of  the  model  by  the  end  of  the  radiative  growth  phase. 

4.  Yield  Rate  Function;  f <t) .  The  normalized  yield 
delivery  rate  function,  f(t),  may  be  any  specified 
(non-negative)  function  of  time.  The  examples  used  here 


assuM  square  pulses  or  triangular  pulses  for  conveniencep 


but  the  model  can  utilize  whatever  time  dependence  is 
desiredp  such  as  the  double— exponential  pulse  used  by 
Pomraning  in  reference  2. 

5.  Numerical  Val ues  of  Parameters.  The  values  of 
parametersp  such  as  air  density  and  opacityp  are  taken  from 
reference  2.  This  assures  comparability  of  the  results  of 
the  two  models. 

!!!9  Motion  of  Material  Inside  the  Fireball .  Heat 
transfer  by  conductionp  convectionp  or  bulk  motion  is 
neglected.  This  is  reasonable  because  of  the  short  time 
scale  and  because  the  fireball  expands  faster  than  the  speed 
of  sound  in  the  fireball  during  the  radiative  phasep  by 
definition. 

Assumptions  Revi sed  from  the  Early  Time  Model 

As  described  in  section  II. Ap  the  major  assumptions  of 
a  model  concern  burnout p  energy  transport  mechanismsp  and 
equilibrium.  The  major  assumptions  needed  for  the  homo¬ 
geneous  models  developed  here  are  presented  in  that  order. 

Burnout .  Difficulties  with  the  early  time  model 
arose  from  the  assumption  that  burnout  is  maintained  until 
the  streaming  radiation  has  been  fully  absorbed.  For  the 
homogeneous  modelsp  it  is  instead  assumed  that  burnout  is 
maintained  only  until  the  absorption  rate  required  to 
maintain  burnout  exceeds  the  ability  of  the  airp  at  burnout 
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opacity,  to  absorb  anergy  from  tha  fiald  of  straaaing 
radiation  prasant.  Tha  datails  of  mathaaatically  aodaling 
this  condition  for  tha  and  of  tha  burnout  subphasa  ara 
davalopad  in  saction  IV. E. 3.  The  assumptions  about  opacity 
of  tha  air  naeded  for  that  daval opment  ara  discussed  next. 

2.  of  the  Air  to  the  Thermal.  Radiation. 

Tha  air  opacity  to  thermal  radiation,  /C  Assumed  to  ba 

an  invarsa  cubic  function  of  tha  energy  of  tha  tharmal 
photons:  Since  the  energy  of  tha  tharmal 

photons  is  proportional  to  the  matin^ial  temperature  of  tha 
air  that  emitted  them,  it  is  equivalently  assumed  that  tha 
opacity  is  an  invarsa  cubic  function  of  the  air  temperature, 
for  tamper aturas  up  to  burnout  temperature.  However,  since 
the  main  interactions  ara  bound-free,  tha  air  is  assumed  to 
ba  transparent  whan  fully  ionized,  so  that  ^  i%  assumed  to 
ba  zero  above  T^.  Thus,  dropping  the  subscript  "t"  for 
reasons  discussed  in  paragraph  3,  below: 

"b'V* 

For  T_  >  T.  (1) 

M  D 

/C(T^)  -  0 

where  means  fC 
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Opacity  of  the  Air  to  Streaming  Radiation. 

As  for  the  thermal  radiation,  the  opacity  to  streaming 
radiation  is  an  inverse  cubic  function  of  the  streaming 
radiation  photon  energy,  hence  of  the  streaming  radiation 
source  temperature,  T  .  However,  T  is  a  complex  function  of 
time  and  is  dependent  upon  the  design  and  construction  of 
the  nuclear  explosive,  and  is  not  presumed  to  be  known.  It 
is  known  that  T  >  T  (otherwise  the  heat  transfer  would 

S  HI 

be  inward).  Therefore,  the  opacity  to  streaming  radiation  is 
less  than  the  opacity  to  thermal  radiation.  The  three- 
temperature  model  uses  this  bounding  value,  as  the 
value  for  /l^^,  which  must  be  recognized  as  an  approximation. 
This  is  a  better  approximation  than  the  early  time  model’s 
assumption  of  infinite  opacity  to  streaming  radiation.  Both 
opacities,  and  will  subsequently  be  referred  to  as 
K  without  subscript  and  without  distinction  between  them. 

4.  Energy  Transport  Mechanisms.  Predictive 
difficulties  with  the  early  time  model  arose  from  the 
assumption  that  a  period  of  streaming  only  transport  is 
followed  by  a  period  of  diffusion  only  transport,  with  no 
period  of  overlap.  Pedagogical  difficulties  arose  from  the 
techniques  used  to  model  these  two  transport  mechanisms, 
namely  the  integrations  over  retarded  time,  and  the  moments 
ecfuations.  To  avoid  these  difficulties,  three  major 
assumptions  are  made. 
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Str»««ii  ng  Mach  an  iM 


It  is  asBuaad  that  tha 


rata  of  transfar  of  anargy  by  atraaalng  through  tha  total 
surfaca  araa  of  tha  firaball  ia  givan  by 

fStraaaing  Rata  through  Surfacal  ■  c  A  (2) 


Mhara,  by  dafinition, 

A  ■  faraa  of  tha  firaball  aurfacal 

■  favaraga  valua  of  tha  straaming  radiation 
dansity  in  tha  firaball! 


(3) 

(4) 


Thia  ia  tha  aaaa  modal  uaad  in  tha  tha  burnout  phaaa  of  tha 
aarly  tima  modal  axcapt  for  tha  additional  aaaumption  that 
E^  ia  uniformly  distributad  throughout  tha  firaball,  ao  tha 
valua  of  E^  at  tha  aurfaca  ia  Juat  tha  avaraga  valua  of  E^ 
inaida  tha  firaball.  Thia  aaaumption  aliminataa  tha  naad  for 
intagrala  ovar  rat ar dad  tima. 

Dif fuaion  Machaniam.  It  ia  aaaumad  that  tha 
corraaponding  anargy  tranafar  rata  for  dif fuaion  of  tha 
tharmal  radiation  ia  givan  by: 

fStraaming  Rata  through  Surfaca!  ■  ^  E^  A  (S) 

Thia  aaaumaa  that  radiation  atraama  at  apaad  c  but  diffuaaa 
at  apaad  c/4  and  ia  aquivalant  to  aaauming  tha  radiatad 
powar  aquala  that  of  a  black  body  in  aquilibrium  at  a 
tamparatura  of  T^^^  aa  ia  ahown  in  aaction  ZV.F. l.b. 
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Overlap  Period.  It  is  assumed  that  streaming 
transport  occurs  as  long  as  there  is  streaming  radiation 
present  in  non-negligible  amounts.  It  is  assumed  that 
di'f  fusion  cannot  transport  energy  through  the  fireball 
surface  to  contribute  to  the  growth  rate  as  long  as  the 
fireball  is  expanding  at  a  rate  greater  than  c/4  due  to 
the  streaming  radiation.  Conceptually,  this  is  because  the 
thermal  radiation  can't  diffuse  outward  fast  enough  to  keep 
up  with  the  surface,  which  is  ever  receding  away  from  it. 
Thus,  there  is  an  overlap  period  of  transport  by  streaming 
and  diffusion  both.  This  overlap  period  starts  whcm  dR/dt 
falls  to  c/4  and  lasts  until  the  streaming  radiation  becomes 
negligible  (in  the  equilibrium  model,  discussed  below). 
Because  of  this  overlap  period,  the  three-temperature 
homogeneous  model  has  smooth  transitions  and  behaves  as  a 
single,  continuous  model,  unlike  the  early  time  model  with 
its  two  distinct  phases. 

S.  Eoui librium.  While  the  early  time  model  never 
explicitly  assumes  equilibrium,  it  will  be  necessary  to  do 
so  for  the  one-temperature  homogeneous  model .  The  approach 
used  is  to  develop  the  three-temperature  homogeneous  model 
from  the  assumptions  presented  above.  The  model  predicts, 
numerically,  that  equilibrium  is  closely  approached  some 
tens  of  shakes  after  the  streaming  radiation  becomes 
negligible.  However,  the  three-temperature  finite 
differencing  scheme  becomes  unstable  when  equilibrium  is 


I 


approachad  too  closely.  (This  instability  is  discussed  in 
the  final  paragraph  of  Appendix  B. )  Therefore,  it  is 
necessary  to  assume  equilibrium  at  this  point,  and  a  one- 
temperature  model  is  used  for  the  remainder  of  the 
computations.  (In  practice,  the  cross-over  point  between  the 
two  numerical  methods  was  established  manually  by 
examination  of  the  output  of  the  finite  differencing 
scheme.)  The  one-temperature  model  incorporates  all  the  same 
assumptions  and  modeling  of  the  physical  processes  that  are 
developed  for  the  three-temperature  model,  so  the  difference 
is  only  of  computational,  rather  than  conce^ptual , 
significance.  Therefore,  the  combination  of  the  two  will  be 
referred  to  as  the  homogeneous  model. 
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IV.  Nonequi librium  Phase; 
Three-Temperature  Homogeneous  Model 


A.  The  Three  Energy  Densities 


The  energy  o4  the  fireball  is  partitioned  among  three 
distinct  forms:  streaming  radiation,  material  heat,  and 
thermal  radiation.  The  three  temperatures  are  T^,  T^  and 
T^py  and  are  defined  next. 


1.  Streaming  Radiation.  The  streaming  radiation  comes 
directly  from  the  nuclear  explosive  device  with  an  energy 
distribution  which  may  or  may  not  be  Plankian.  Since  this 
analysis  uses  the  (bounding)  value  of  the  fireball  opacity 
to  thermal  radiation  as  an  approximation  for  the  opacity  to 
streaming  radiation  (see  section  III.C.3),  the  source 
temperature  of  the  streaming  radiation,  T^,  is  not  of 
concern.  Therefore,  the  streaming  radiation  will  be 
characterized  only  by  its  energy  density,  E^. 

2.  Material  Energy.  The  material  energy  is 
characterized  by  its  thermodynamic  temperature,  T  .  The  ions 

III 

and  electrons  of  the  fireball  plasma  are  assumed  to  be  in 
thermodynamic  equilibrium  (with  each  other),  so  distinct  ion 
and  electron  temperatures  are  not  considered.  Since  the  heat 
capacity  is  assumed  constant,  the  energy  density  of  the 
material  is: 


E 


m 


'=vTm 


(6) 
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3.  Thermal  Radiation.  The  thermal  radiation  will  not 
be  in  black  body  equilibrium  during  this  phase,  but  is 
approaching  such  equilibrium.  Therefore,  it  is  reasonable  to 
characterize  the  thermal  radiation  field  by  the  temperature 
of  a  Plankian  distribution  of  equal  energy  density.  This 
"effective  temperature"  is  defined  by  the  relation: 


B.  Physical  Processes 

The  physical  processes  considered  in  the  model  are 

summarized  here  and  detailed  in  the  rest  of  this  section. 

1.  Streaming  B§dlation  Processes.  The  streaming 

radiation  is  emitted  from  the  nuclear  explosive  device.  It 

enters  the  fireball  at  its  interior  boundary.  It  is  subject 

to  absorption  within  the  fireball  as  it  streams  radially  to 

the  exterior  boundary.  At  the  surface  of  the  fireball,  the 

streaming  radiation  enters  the  (cold)  ambient  air,  and  is 

absorbed  in  a  short  distance,  on  the  order  of  a  mean  free 

path  of  an  x-ray  in  ambient  air  (a  few  cm).  Thus,  one  may 

visualize  a  thin  absorption  layer  on  the  surface  of  the 

fireball.  The  ambient  air  ' n  this  layer  is  heated  to  the 

fireball  material  temperature,  T  ,  ionizing  it  and  advancing 

ni 

the  surface  of  the  fireball  outward.  In  a  time  increment  dt, 

the  radius  of  the  fireball  is  increased  by  dR,  producing  an 

additional  volume,  dV,  of  air  heated  to  T  and  filled  with 

in 

streaming  energy  at  density  E^.  The  shell  of  thickness  dR 
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Nhich  constitutes  dV  is  not  the  thin  absorption  layer,  but 
rather  is  the  volume  swept  out  by  the  advancing  surface  of 
the  fireball  in  the  time  increment  dt.  Since  the  eajor 
assumption  of  this  model  is  that  temperatures  are 
homogeneous  inside  the  fireball  and  zero  outside  the 
fireball,  this  requires  the  idealization  that  the  absorption 
layer  is  of  zero  thickness.  In  recognition  of  this,  this 
absorption  of  radiation  will  be  referred  to  as  absorption  at 
the  surface  of  the  fireball. 

2.  Material  Energy  Interactions.  Because  the  period 
modeled  here  is  entirely  before  hydrodynamic  separation, 
heat  conduction,  convection,  or  other  direct  interaction  are 
neglectedi  radiative  interactions  only  are  considered.  These 
interactions  occur  both  within  and  at  the  advancing  surface 
of  the  fireball.  Within  the  fireball,  material  energy  is 
produced  by  absorption  of  streaming  and  thermal  radiations! 
material  energy  is  lost  by  emission  of  thermal  radiation. 
These  interior  effects  tend  to  raise  or  lower  the 
temperature,  T  .  As  the  surface  advances,  material  energy  is 
produced  by  absorption  of  streaming  and  thermal  radiations, 
causing  fireball  growth  rather  than  temperature  changes. 

Thermal  Processes.  The  thermal  radiation 
is  formed  by  emission  throughout  the  volume  of  the  fireball. 
It  is  also  absorbed  within  that  region.  Thermal  radiation 
moves  to  the  surface  of  the  fireball  by  diffusion  and  is 
absorbed  there. 
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C.  Total  Enorqy  Bal ance  Equations 

1.  Bansral  B1 obal  Bal ance  Equation.  Tho  dlf forontlal 
aquations  of  the  three-temperature  model  will  be  developed 
by  means  of  global  energy  balances  on  each  of  the  three 
energies.  By  global  balance,  1  mean  the  energy  density 
Integrated  over  all  space  (the  total  energy,  Q)  balanced 
before  and  after  a  time  Increment,  dt.  Such  balances  are 
simplified  by  the  assumption  that  the  energy  densities  are 
spatially  uniform  Inside  the  fireball  and  zero  outside  the 
fireball.  Thus: 

E  {Total  Energy  of  form  x  at  time  t> 

dV  -  E,(t)  (8) 

where  "x"  Indicates  streaming  or  material  or  thermal. 

The  global  energy  balance  equations  are  of  the  form 

"New  Value  *  Old  Value  Production  In  dt  -  Loss  In  dt" 

It  Is  useful  to  distinguish  between  product Ion /loss 
occurring  In  the  fireball  Interior  and  product Ion /loss 
occurring  at  the  surface  because  of  the  different  nature  of 
these  processes.  With  this  distinction,  the  general  global 
balance  equation  becomes: 
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{Total  anargy  at  tiaa  t-*-dt>  ■  (9) 

{Total  anargy  at  tima  t> 

{Enargy  producad  insida  during  dt> 

-  {Enargy  abaorbed  insida  during  dt> 

{Enargy  produced  at  surface  during  dt> 

-  {Energy  absorbed  at  surface  during  dt} 

Mhera  *m*  indicates  streaming  or  material  or  thermal. 

2.  Production  and  Loss  at  Surface  of  Firebal 1 .  By 
definition,  during  time  interval  dt,  the  fireball  grows  in 
volume  by  the  addition  of  a  shell  of  thickness  dR  and  volume 
dV  onto  the  original  fireball  of  radius  R  and  volume  V, 
producing  a  new  radius  R*  and  volume  V*  at  time  t*  ■  t-*-dt. 
The  conversion  of  energy  at  the  surface  during  dt  can  most 
easily  be  visualized  using  these  terms,  and  is  described 
next. 

For  streaming  and  thermal  radiation  energies, 
there  is  no  production  at  the  surface.  Absorption  occurs  on 
the  surface,  which  moves  from  R  to  R'  during  dt.  This 
absorption  is  found  by: 
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{Radiation  anergy  absorbed  at  sur-face  during  dt)  *■ 
{Radiation  energy  Nhich  enters  dV 

through  sur-face  at  R  during  dt>  (10) 

-  {Radiation  energy  which  is  still  present 

in  dV  as  radiation  energy  o-f  type  r  at  time  t'> 

where  *r*  indicates  streaming  or  thermal. 

For  material  energy,  there  is  no  loss  on  the  sur-face. 
The  production  on  the  sur-face  (by  absorption  o-f  streaming 
and  thermal  radiation)  is  found  by: 

{Material  energy  produced  at  surface  during  dt>  > 

{Material  energy  Q  present  in  dV  at  time  t'>  (11) 

fll 

3.  Production  and  Loss  Inside  the  Pi  rebel,!.  The 
average  volume  of  the  fireball  during  dt  is: 

-  V  +  dV/2  (12) 

and  the  interior  conversion  of  energy  can  be  visualized  as: 

^Production  energy  Qx  inside  fireball  during  dt>  ■ 
^rSdCJtiSI  -n-rQy  density  Em  at  t>  *  *  dt 

where  *)<*  indicates  streaming  or  material  or  thermal,  and 
the  absorption/production  rate  is  assumed  constant  during 
dt. 
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Next,  th«  various  production  and  loss  aschanisMs  Hill 
bs  considorsd  and  rats  squat ions  dovslopsd.  Tho  rssults  will 
than  bo  substitutsd  in  tho  global  snorgy  balancss  to  obtain 
tho  dosirod  difforential  ocpiations  and  finito  di-fforonco 
aquations. 

D.  Cantor  Point  Boundary  Conditions 

Tho  anargy  rolaaso  is  assunod  to  ba  ^roo  a  dovica  of 

nogligibla  volumo  (a  “point  sourca").  Thorofora,  tho 

intarior  boundary  raducas  to  a  cantor  point  for  tha  purposa 

of  spacifying  boundary  conditions.  These  boundary  conditions 

are  simplified  by  the  assumptions  that  each  energy  density 

is  spatially  uniform  within  the  fireball  (S-temperature- 

homogeneous)  and  that  the  material  and  the  thermal  radiation 

do  not  interact  directly  with  tha  axplosiva  device. 

Consaquently,  E  and  automatically  meet  the  requirements 

of  finitenass,  continuity,  and  no  net  flux  through  tha 

interior  boundary.  Thus,  no  intarior  boundary  conditions 

need  ba  explicitly  specified  few  E  and  E.  . 

m  X, 

Tha  streaming  radiation  enters  tha  fireball  from  the 
point  sourca.  Tha  associated  boundary  condition  can  ba 
treated  as  a  term  in  tha  global  balance  equation  for  tha 
streaming  radiation  anargy  as  follows: 

{Streaming  anargy  produced  inside 

-  Y  f(t)  dt  (14) 

(at  canter  of)  fireball  during  dt> 
where  Y  is  tha  Yield. 
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E.  Eouili brut ion  Relations  Within  the  Fireball 

E«i eei on  of  Thermal  Radiation.  The  thermal  radia¬ 
tion  is  emitted  at  a  rata  proportional  to  the  density  of 
air,  p,  its  opacity,  K*  (which,  as  described  in  section 
XII.C,  is  assumed  to  be  an  inverse  cubic  function  of  its 
teoperatur  e)  and  to  the  temperature  to  the  fourth  power  as 
follows: 


5t  ^ 


emi ssion 


-  C(T,)  aT 

In  11 


where 


(IS) 


so  that 


C(T,)  -  pcAC^  <VT.>’ 


(16) 


{Thermal  radiation  energy  produced  in  fireball  during  dt> 

-  <=<%>  dt  <i7) 

■  {Material  energy  lost  inside  fireball  during  dt> 


2.  Absorption  of  Ihermal  Radiation.  The  thermal  radia¬ 
tion  is  absorbed  at  a  rate  proportional  to  the  density  of 
thermal  radiation  present  to  be  absorbed,  and  with  the  same 
proportionality  coefficient  as  for  emission.  Thus: 


and 


5t  ^t 


-  -  C(T  )  E, 


absorption 


(18) 


{Thermal  radiation  energy  absorbed 
inside  fireball  during  dt> 


C(TJ  E^ 


Ve  *** 


(19) 
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3>  Ab»carption  of  Streaming  Radiation.  A  portion  of  tho 
■traaMing  radiation  is  absorbed  by  the  ionized  air  in  the 
fireball.  The  rate  of  absorption  depends  on  the  conditions 
in  the  fireball  in  the  following  way: 

If  T  -  T.  and  CaT^  -  E.  3  <  then 

fli  D  O  ^  V 

fStreaaing  energy  absorbed  inside  fireball  during  dt> 

-  dt 

otherwi se 

<6treamlng  energy  absorbed  inside  fireball  during  dt> 

-  ^  ''.v.  «  '2®> 

In  general  terms,  there  are  two  criteria  to  consider: 
the  streaming  radiation  can  be  intense  or  weaki  the  fireball 
can  be  at  burnout  or  below  burnout.  If  the  fireball  is  below 
burnout  and  the  streaming  radiation  is  weak,  then  the 
absorption  of  streaming  radiation  will  be  less  than  the  net 
loss  of  E  by  emission  of  thermal  radiation.  The  fireball 
will  be  cooling. 

Conversely,  if  the  fireball  is  below  burnout  and  the 
streaming  radiation  is  intense,  then  the  absorption  of 
streaming  radiation  will  exceed  the  net  loss  of  E^  by 
emission  of  thermal  radiation.  T_  will  be  rising. 

ill 

But  if  the  fireball  is  already  at  burnout,  then  it 
cannot  absorb  more  streaming  radiation  than  it  loses  to 
thermal  radiation,  because  that  would  raise  its  temperature 


•bow  burnout,  which  would  «•!<•  it  tr*nspAr«nt,  proventing 
absorption  o-f  streaming  radiation.  So  no  matter  how  intense 
the  streaming  radiation,  the  fireball  can  only  absorb  enough 
to  maintain  T^,  i.e.,  enough  to  equal  the  net  loss  of  by 
emission  of  thermal  radiation. 

The  final  case  to  consider  is  if  the  fireball  is  at 
burnout  and  the  streaming  radiation  is  weak.  In  this  case, 
the  absorption  of  streaming  radiation  is  small  and  the 
fireball  cools,  as  in  the  first  case. 

Zn  order  to  quantify  these  arguments,  consider  first  the 

cases  where  T  <  T.  .  Since  the  same  value  of  opacity  is 

in  D 

being  assumed  for  the  streaming  radiation  as  for  the  thermal 
radiation,  the  absorption  rate  of  streaming  radiation 
density,  E^,  depends  upon  the  air  density,  the  opacity,  and 
upon  E^  in  the  same  way  as  shown  above  in  the  case  of 
absorption  of  thermal  radiation.  Thus,  the  rate  of 
absorption  of  streaming  radiation  is  given  by: 


{Total  Streaming  Energy,  Q^, 
absorbed  inside  fireball  during  dt> 


^  Ve 


(21) 


Zn  the  case  where  T  «  T.  and  the  streaming  radiation  is 
weak,  this  same  relation  holds.  Zf  the  streaming  radiation 
is  intense,  however,  the  absorption  of  streaming  radiation 
cannot  be  large  enough  to  cause  the  air  temperature  to  rise, 
as  described  above.  Thus,  the  rate  of  absorption  of 
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•tr»a«ing  radiation  cannot  exceed  the  amount  needed  to 
maintain  the  air  temperature  at  burnout,  which  is  equal  to 
the  net  loss  rate  of  material  energy  (by  emission  and 
absorption  of  thermal  radiation).  Hence,  the  absorption  rate 
of  streaming  radiation  is  bounded  by: 


CNet  loss  of  Q  by  emission  of  . 

®  -  C(T  >  Car;  -  E  3  V  dt 

inside  fireball  during  dt>  m  m  c  ave  ^^2) 


Thus, 

then 

the 

absorption  is  given 

by 

equation 

(21) 

unless  it 

would 

excc»ed  the  bound  givcm 

by 

equation 

(22), 

i.e.  unless 

CaT^- 

m 

<  E.. 

Combining  equations  (21)  and  (22),  together  with  their 
conditions  of  applicability,  gives  the  result  initially 
stated,  equation  (20). 


4.  Production  of  Material  Energy.  Material  enewgy  is 
produced  inside  the  fireball  by  absorption  of  both  streaming 
and  thermal  radiation.  Combining  equations  (19)  and  (20) 
gives  the  term  needed  for  the  global  balance  equation  for 
material  energy: 


If  T  -  T.  and  CaT^  -  E.  3  <  E  then 
mo  m  t  s 

{Material  energy  produced  inside  fireball  during  dt> 

■  {Material  energy  absorbed  inside  fireball  during  dt> 
other wi se 

{Material  energy  produced  inside  fireball  during  dt> 


C(T  ) 

in 


CE. 


Ve 


(23) 
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F.  Firgball  Growth :  Surface  Boundary  Conditions 

1.  Energy  Transport.  Energy  is  transported  to  the 
surface  of  the  fireball  by  two  mechanises  of  radiatiyo 
transport,  streaming  and  diffusion. 

Streaming.  Transport  of  the  streaming 

radiation  causes  energy  to  arriye  at  the  surface,  radius  R, 
at  a  rata  of  cE^  per  unit  area,  since  E^  is  assumed 
uniform  throughout  the  fireball.  Thus,  in  terms  of  a  global 
energy  balance,  the  total  rate  of  energy  arriying  at  the 
surface  of  the  fireball  is  cE^A  where  A  is  the  area  of  the 
surface,  and 

{Streaming  radiation  energy 

which  enters  dV  through  •  c  E^  A  dt  (24) 

surface  at  R  in  time  dt> 

Diffusion.  Transport  of  the  thermal  radiation 
is  by  diffusion.  This  process  is  essentially  that  of  a 
radiating  black  body  at  a  temperature  of  T^^ ^ ,  which  would 
radiate  with  power  per  unit  area  of  (J  T^^^  Since 
O'  *  ^  •  the  total  rate  of  arriyal  of  thermal  radiation  at 
the  (stationary)  surface  is 

■  I 


{Thermal  radiation  energy 
which  enters  dV  through 


(26) 


2.  Fireball  Growth.  The  growth  o-f  the  fireball  is 
determined  by  an  energy  balance  on  the  incremental  volume, 
dVy  burned  out  in  time  increment  dt.  There  are  three 
possible  cases  to  consider.  These  consist  of  growth  caused 
by  streaming  only,  by  diffusion  only,  and  by  both  together. 

Growth  by  Streaming  Qn^y.  In  this  case,  the 
energy  balance  requires  that  the  energy  entering  dV,  by 
streaming  through  the  surface,  must  equal  th'e  energy  density 
in  dV  times  dV.  Thus,  since  only  streaming  radiation  enters 
dV: 


Solving: 


c  E  A  dt  »  <E  +  E_)  dV 

s  s  m 


(27) 


•  A  c  E_  /  (E_  ♦  E.)  (28) 

CIb  B  B  fli 


This  case  applies  when  the  streaming  radiation  is  so  intense 
as  to  cause  dR/dt  >  c/4  so  that  the  diffusion  of  thermal 
radiation  essentially  cannot  keep  up  with  the  advancing 
surface  of  the  fireball.  This  is  an  approximation  since  some 
thermal  photons  must  be  directed  nearly  directly  outward, 
but  this  approximation  will  prove  to  be  a  useful  one. 

Growth  by  Diffusion  Only-  In  this  case,  the 
energy  balance  is  the  same  as  for  streaming,  but  E^  is 
replaced  by  E^  and  c  by  c/4,  since  radiation  streams  at  c 
but  diffuses  at  c/4.  Thus: 

-  E^  A  dt  -  <E^  +  E^)  dV  (29) 
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Solving: 

dt  ■  4  ^ 

This  c«SB  applies  when  the  streaming  radiation  density 
becomes  negligible,  »iihich  is  somewhat  later  than  the  end  of 
the  yield  pulse  which  produces  the  streaming  radiation. 

c-  §CS!lS^tl  ky  Diffusion  and  Streaming  Superposed. 
In  this  case,  the  energy  balance  on  dV  assumes  that  dV  is 
filled  with  streaming  and  thermal  radiations  at  the  same 
densities  as  in  the  interior  of  the  fireball.  The  excess  of 
energy  entering  dV  is  present  as  material  energy  at  the  same 
temperature  as  inside  the  fireball.  Thus: 

A  <c  +  I  E^)  dt  •  (E^  +  E^  ♦  E^)  dV  (31) 

Solving: 

“  -  A  c  <E.  +  E. /4)  /  (E_  +  E.  +  E_)  (32) 

This  case  applies  when  dR/dt  <  c/4  and  both  streaming  and 
thermal  radiations  are  present.  Note  that  if  dR/dt  ■  c/4, 

then  both  equation  (2S)  and  equation  (32)  give  the  sanv 
result,  even  if  E^  >0.  Thus,  the  change  from  streaming  only 
to  both  streaming  and  diffusion  is  a  smooth  one.  Note  also 
that  equations  (29)  and  (30)  are  special  cases  of  equations 
(31)  and  (32)  with  E^  >  0. 
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B.  The  Three-Temperature  Mode^  System  of  Equations 

The  model  leads  to  a  system  o-f  six  equations  in  six 
unknonnns:  E^,  E^y  E.  *  V,  A,  and  R.  The  first  equation  is 
obtained  from  a  total  energy  balance  on  dV,  and  is  the 
equation  for  dV/dt.  This  is  really  the  combination  of 
equations  <28),  (30)  and  (32)  as  described  in  section  F, 

above.  Three  more  equations  are  obtained  from  the  global 
energy  balances  described  in  section  C,  above.  The  remaining 
equations  are  the  geometric  relationships  among  volume, 
radius  and  area.  These  geometric  relationships  are 
summarized  first,  then  the  balance  equations  will  be 
obtained.  The  system  of  balance  equations  is  summarized  in 
Table  1)  the  system  of  differential  equations  is  summarized 
in  Table  2. 

Geometric  Equations.  For  a  free  air  burst,  the 
fireball  is  spherical,  so  that: 

\ 

A-47fR^  V-|7fR^-AR/3  (33) 

Zt  is  North  noting  at  this  point  that  none  of  the 
derivation  of  this  model  has  been  dependent  on  the 
assumption  of  spherical  geometry.  Therefore,  if  a  ground 
burst  were  to  be  analyzed,  it  would  only  be  necessary  to 
assume  that  no  energy  is  transferred  into  the  ground  and  to 
replace  equations  (33)  with  those  appropriate  to  the 
geometry  considered. 


2.  Streawinq  6sdi.atign  Energy  Balance.  The  global 
balance  equation  is: 

{Streaming  Energy  at  t*>  (74) 

{Streaming  Energy  at  t> 

{Production  from  Yield  during  dt> 

-  {Absorption  inside  fireball  during  dt> 

-  {Transfer  through  surface  at  R  during  dt> 
{Streaming  Energy  present  in  dV  at  t'> 


From  equation  <8) 

{Streaming  Energy  at  t>  ■  V  E^(t) 
and  similarly 

{Streaming  Energy  at  t’>  «  V  E  (t*> 

s 

Recalling  equation  (14) 

{Production  from  Yield  during  dt>  ■  Y  f(t)  dt 
Defining  the  term  dE  to  be  the  energy  density  converted 

ShI 

from  streaming  to  material  energy  density  during  dt: 

{Absorption  inside  Fireball  during  dt>  ■  '^ave  **^sm 
Recalling  equation  (24) 

{Transfer  through  surface  at  R  during  dt>  ••  c  E^  A  dt 
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Assumi ng 

the 

stream! ng 

energy 

density 

in 

dV  at 

t' 

is  the 

same  as 

the 

streaming 

energy 

density 

in 

V  at 

t 

(as  was 

tacitly  assumed  in  the  derivation  o-f  dV/dt),  then 


(Streaming  Energy  present  in  dV  at  t'>  ■  dV 

Substitution  o-f  these  expressions  into  equation  (34)  gives 
the  desired  balance  equation: 


E^»  V'  «  E^  V  +  Y  •f<t)  dt 

-  c  E  A  dt 

B 


ave 


E  dV 

s 


(35) 


'  “T.  and 

CaT^-E. 3 

< 

E 

m  b 

m  t 

B 

Then 

dE 

m 

C(T_) 

CaT^-E.!  dt 

sm 

m 

m  t 

OtherMi se 

dE 

sm 

- 

C(T  > 

In 

E^  dt 

(36) 

Then,  since 

ZZ  <E  V)  -  limit  EE  <t+dt) V(t+dt)  -  E, (t) V(t> 3/dt  (37) 
•  dt>0  •  * 


Equation  (3S)  implies: 


d 

dt 

(E^V) 

-  V 

dE  /dt  + 

Sin 

* 

dt 

Y  f(t)  -  c  E^  A 

(38) 

where 

If 

and 

CaT^~E.3  < 

HI  V 

Then 

- 

C(T.) 

in 

Otherwi se 

dEsm/dt 

fli 

C(T_) 

n 

(39) 
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Material  Energy  Balance.  The  global  balance 

equation  is: 

{Material  Energy  at  t'>  ■  ^4^, 

{Material  Energy  at  t> 

♦  {Gain  by  absorption  of  radiation  inside  during  dt> 
{Loss  by  emission  of  radiation  inside  during  dt> 

♦  {Material  Energy  in  dV  at  t'> 

Defining  the  term  dE^^  to  be  the  net  energy  density 
converted  from  material  energy  density  to  thermal  radiation 
density  during  dt: 


{Gain  by  absorption  of  thermal  radiation  inside  during  dt> 
”  {Loss  by  emission  inside  during  dt>  (41) 

“  -  dE  .  V 

mt  ave 

Adding  the  gain  by  absorption  of  streaming  radiation: 


{Gain  by  absorption  inside  during  dt> 

-  {Loss  by  emission  inside  during  dt> 


(42) 


-  ''.V.  -  «.t> 


The  other  terms  parallel  those  for  streaming  radiation,  so 
the  balance  equation  is: 


E.'  V'  -  E.  V 


avi 


•  "^mm  -  ♦  E.  dV 


(43) 
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Mhsroi 


«  ««> 

Nhich  leadm  to  tho  di-f-farontial  equation: 

<VE^>  -  V  CdE__/dt  -  dE  . /dtD  +  E  <^S> 

flu  fll  SA  Aw  n  Clw 


Nhere: 

dE^/dt  -  C<T^>  C.t;-E^3  C46) 

4.  Thermal  Radiation  ^ergy  Balance.  The  global 
balance  equation  is: 

<Thermal  Radiation  Energy  at  t'>  <■  (47> 

{Thermal  Radiation  Energy  at  t> 


*  {Production  by  emission  inside  during  dt> 

-  {Loss  by  absorption  inside  during  dt> 

-  {Transfer  through  surface  at  R  during  dt> 

{Thermal  Radiation  Energy  present  in  dV  at  t'> 

As  Mas  discussed  in  section  F*  above*  there  are  two  cases  to 
consider  here.  If  dR/dt  >  c/4*  then  the  thermal  radiation 
cannot  enter  dV  in  sufficient  quantity  to  fill  dV  with 
thermal  radiation  of  density  E^.  It  is  therefore  assumed 
that  no  conversion  of  thermal  radiation  to  material  energy 
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occurs  in  this  cass.  As  a  rssult. 

If  dR/dt  >  c/4  than 

CTransfsr  through  surfacs  at  R  during  dt>  ■ 

(48) 

fThsrmal  Radiation  Ensrgy  prsssnt  in  dV  at  t*> 

OthSTMiss 

{Transfsr  through  surfacs  at  R  during  dt>  »  ^  ^4.  A 

^  '  (49) 

fThsraal  Radiation  Ensrgy  prsssnt  in  dV  at  t'}  ■  dV 

Substituting  aquations  (8),  (41),  (48)  and  (49)  into  (47) 

givss  ths  global  balancs  aquation  for  tharsal  radiation: 

If  dR/dt  >  c/4 

Than  E^'  V'  -  E^  V  ^  dE^^  (50) 

Othsrwi sa 

E^'  V'  -  E^  V  *  dE^^  -  I  E^  A  dt  *  E^  dV  (51) 

Taking  tha  limit  as  dtM),  and  substituting  aquation  (44)  for 
dE^^/dt  givss  the  differential  aquations  for  the  tNO  cases: 

If  dR/dt  >  c/4 

Than  (VE. )  -  V  C(T_)  CaT^  -  E.  3  (52) 

Ov  b  HI  01  H 

OthsTHi sa 

3t  "^t*  -  ''  '*^2  -  ♦  =t  5?  -  I  S 
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A  numarical  approach  to  solving  thsss  squat! on  Mill  bs 
proBsntsd  a-ftor  considoring  ths  initial  and  final  conditions 
■for  ths  solution. 

H.  Initial,  Conditions 

Ths  above  system  of  equations  is  to  be  solved  with  the 
following  initial  conditions: 


(94) 


Because  the  streaming  radiation  starts  at  a  point,  it  has  an 
initially  infinite  value,  causing  dR/dt  ■  c. 


I-  Final  Conditions 

The  three-temperature  model  applies  until  two 
conditions  are  both  met:  the  streaming  radiation  has  become 
negligible,  and  the  material  and  thermal  radiation  energies 
are  in  equilibrium,  i.e. : 


and 


(99) 


Once  these  conditions  are  met,  the  one-temperature  model 
developed  in  the  next  section  may  be  used. 
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TABLE  1 


XbE5S  Temperature  Model  Balance  Equation* 


A  -  4  7f  R" 


V  -  ItTr’  -  AR/3 


E  •  V'  -  E.  V  +  Y  fit)  dt  -  dE  V 
•  •  Bffl  ave 

-  c  E_  A  dt  +  E  dV 


E  '  V*  - 


V  +  <‘1E^^  -  dE  .  )  +  E  dV 

m  ave  em  mt  m 


I-f  dR/dt  >  c/4  then 


E^'  V* 


E.  V  +  dE  .  V 
t  At  ave 


and 


c  E.  A  dt  -  <E.  +  E  )  dV 

m  «  m 


other Hi se 


and 


^t'  V'  -  E,  V  ^  dE„,  .  I  E^  A  dt  ^  E^  dV 


A  <c  E.  +  ;  E.  )  dt  »  <E_  +  E.  +  E  )  dV 


Mhere: 


dE,t  -  C(T.)  C.T;-E^3  dt 


and 


■  “T.  and 

CaT  -E^l  < 

E 

A  b 

At 

s 

Than 

<*E—  ■ 

C(T  ) 

CaT^H 

SA 

A 

A 

Otherni se 

C(T_) 

in 

(33) 

(39) 

(43) 

(50) 

(27) 

(51) 

(31) 

(44) 

(36) 
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TABLE  2 


Thr«»  Taapcwature  Model  Differential  Equatione 


A  -  47f  V  -  ItT  R*  -  A 

R  /  3 

(33) 

dt 

-  V  dE_  /dt  +  S  +  Y  f (t)  - 

-  c  E,  A 

(38) 

5t 

-  V  CdE  /dt  -  dE^^/dtl  +  E  ~ 

Ml  Au  m  Ow 

(45) 

If  dR/dt 

>  c/4  then 

gj  (VE^>  -  V  C<T„)  I.t;;  -  E^l 

(52) 

and 

K  -  A  =  E,  /  <E.  -  E.> 

(28) 

Otherwise 

at  ‘'^t*  ■  ''  ‘='T.'  ^ 

dV  C  e- 
dt  “  4  H 

A 

(53) 

and 

~  -  A  c  <E_  +  E. /4)  /  (E_  +  E*. 

(32) 

Where: 

dE  ./dt  »  C<T_)  CaTf“E*3 

fli  w  IS  W 

(46) 

and 

If  T  -T.  and  CaT^-E.l  <  & 

AD  A  V  • 

Then  dE  /dt  ■  C(T^) 

Ml  A 

(39) 

otherwise  dE _ /dt  ■  C(T^) 

AA  A 

^s 

J.  Finite  Differencing  the  Three-Temperature  Model 

This  section  describes  an  algorithm  for  solving  the 
three-temperature  model  equations  by  finite  differencing.  A 
program  implementing  this  method  in  Microsoft  Basic  version 
5.2  is  listed  in  Appendix  B. 

1 •  Evaluate  dV. 

a.  Assuming  streaming  only,  solve  equation  (27): 

dV  -  c  A  dt  E  /  (E_  +  E_)  (56) 

m  mm 

b.  Verify  that  dR/dt  >  c/4  using: 

dR  «  3  dV  /  A  (57) 


c.  If  dR  <  (c/4)dt  then  streaming  and  diffusion 
both  contribute  to  fireball  groMth.  Solve  equation  (31)  for 
dV: 

dV  -  c  A  dt  (E  +  E./4)  /  (E_  +  E.  E^)  (58) 

d.  The  new  volume  is 

V’  ■  V  +  dV  (59) 
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2m  Evaluate  Energy  Penalty  Converaiona.  The  converaion 
aoKing  energy  denaitiea  la  determined  by  equationa  (36)  and 
(44)  aa  followa: 

a.  Find  the  denaity  converted  -from  atreamlng  to 

material  energy,  (^E  ,  aaauming  it  ia  not  bounded  by 

Sin 

tranaparency  o-f  the  air  above  uaing  equation  (36): 

dE.„  -  C(T,)  E.  dt  (60) 

mm  m  S 

b.  Find  the  net  denaity  converted  from  material  to 

thermal  radiation  energy,  uaing  equation  (44): 

*  C(T„)  CaT^I  -  E.  3  dt  (61) 

nz  in  n  ti 


c.  If  the  fireball  ia  at  burnout,  i.e.  if  T. 

HI  D 

then  bound  dE__  by  dE  .  ,  hence  applying  the  conditional 
Sin  mu 

portion  of  equationa  (36).  Thua: 


If 

**E«t  < 

‘‘^am 

Then  let 

‘‘^mt 

Solve  Global  Balance  Equationa  for  Energy  Denaitiea. 
a.  Find  E^’  by  aolving  the  energy  balance, 
equation  (35),  for  E^* : 

“  ‘^^*m  ''ave^  '''  Yf(t)dt/V'  -  c  E^  dt  A/V'  (63) 


SO 


b.  Find  E  *  by  solving  ths  snsrgy  bslancs* 


squstion  (43)  y  -for  E  * : 


E.  ♦  (dE^  -  dE.^)  /  V 


c.  Find  E^'  by  solving  squstion  (50)  or  (51 >y 
dspsnding  on  Mhsthsr  dR/dt  >  c/4  or  not,  rsspsctivsly: 


X«  dR/dt  >  c/4  thsn 


E^  V  /  V  *  dE^^  /  V' 


OthSTMiSS 


V  -  ^  ''.V.  '  ''■  -  <**> 

gvslusts  ths  Rsssining  Vsriablss.  Sol vs  squstions 
(6),  (7),  and  (33)  slgsbrsicslly  for  »  ^sff'* 

in  tsrss  of  E  * ,  E.*,  and  V* ,  rsspecti vsly.  This  complstss 

ths  solution  for  ons  tiss  stsp,  so  all  that  rssains  is  to 


tost  against  ths  stopping  conditions  prsssntsd  in  ssction 
IV. I,  abovs,  and  itorats  as  rsquirsd. 


5.  Initializing  ths  Itsrativs  Procsdurs.  Ths  initial 


conditions  aro  prsssntsd  in  ssction  IV. H, 


abovs.  Ths 


itsrativs  squations  cannot  havs  V  ■  O,  howsvor.  This  is 
avoidsd  by  using  a  first  tiss  stsp  of  sizo  2«dt,  shich, 
togsthsr  with  dR/dt  *  c  givss  a  starting  radius  of: 


R  -  2  c  dt 


wmmm 


Th»  voluns  and  area  corresponding  to  this  radius  are  then 


calculated.  Then,  since  "*  ^  (0)  ■  0 

a  total  energy  balance  is  used  to  -find  E^: 

Y  f(dt)  2  dt  «  V  CE  -i-  E.  1  (68) 

S  D 

These  values  are  then  treated  as  initial  conditions  for  the 
iterative  algorithm. 

Other  Considerations.  The  iterative  scheme 

described  here  implicitly  assumes  that  dV<<V  and  dR<<R. 
If  a  constant  time  step,  dt,  is  used,  then  these  conditions 
are  not  met  for  the  early  time  steps.  Equations  (56),  (57) 
and  (58)  for  dV  and  for  dR  must  be  altered  to  include 
spherical  divergence  in  these  cases  in  order  to  achieve 
accurate  results.  In  some  places,  formulas  must  be 

rearranged  to  avoid  loss  of  precision  from  subtraction  of 
nearly  equal  quantities,  or  to  avoid  fatal  program  errors 
(such  as  division  by  zero).  The  use  of  a  first  time  step  of 
2  dt  rather  than  dt  avoids  such  problems  during  the 
second  time  step,  for  example.  The  program  listed  in 
Appendix  B  incorporates  these  requirements,  and  is  annotated 
with  remark  statements  where  they  occur. 
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Eoui librium  Phase: 
One-Temperature  Homogeneous  Model 


Additional  Assumptions 

The  assumptions  required  -for  the  three-temperature 
model,  detailed  in  section  III  above,  are  also  required  for 
the  one-temperature  model  presented  in  this  section.  In 
addition,  the  one-temperature  model  assumes:* 

—  the  yield  pulse  has  been  completed 

—  the  streaming  radiation  has  been  absorbed 

—  the  thermal  radiation  is  in  equilibrium  with  the 
material  energy. 

Under  these  conditions,  the  energy  densities  inside  the 
fireball  are  characterized  by  a  single  temperature,  T, 


defined  by: 


^m  - 


This  equilibrium  can  never  actually  be  reached,  but  it  is 
approached  closely.  For  example,  50  shakes  after  the  start 
of  a  20  kiloton  burst  of  10  shakes  duration,  the  three- 
temperature  model  predicts  that  T_, ,  and  T_  differ  by  0.3%. 

It  is  worth  noting  at  this  point  that  it  is  not 
possible  to  use  the  computer  implementation  of  the  three- 
temperature  model  beyond  the  point  where  equilibrium  is 
(nearly)  reached,  because  the  numerical  method  becomes 
unstable.  This  difficulty  is  discussed  in  Appendix  B. 
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Th»re<for«,  us«  of  the  one-temperature  model  is  required  in 
order  to  extend  the  computations  to  hydroseparation. 

B.  Approach  to  the  Problem 

The  unknowns  o-f  the  three-temperature  model  are  ^ , 

E. ,  E_y  V,  R,  and  A.  The  approach  to  the  problem 

Ml  s  w  In 

used  here  is  to  express  all  these  unknowns  in  terms  of  a 
single  parameter,  T.  Then  a  differential  equation  for  T(t) 
in  terms  of  physical  constants,  yield,  and  T,  will  be 
obtained.  This  equation  can  be  solved  for  T(t>  using 
conventional  numerical  techniques.  All  the  desired  unknowns 
can  then  be  calculated  from  T. 


C.  Reduction  to  a  Si.ng^e  Parameter!  J 


T^^.y!  from  equation  <69) ! 


%ff  -  T 


T_  -  T 


T  :  from  equation  (69): 


E^:  by  assumption: 


E^:  from  equations  (7)  and  (70): 


^ 


E^:  from  equations  (6)  and  (71): 


^  -  =vT 


(70) 

(71) 


(72) 


(73) 


(74) 


54 


V:  sine*  all  the  yield  has  been  delivered: 

V  (E_  +  E.  +  E  )  -  Y  (75) 

S  b  ID 

Substituting  equations  (72)  through  (74)  into  (75): 

Y 

V  -  — - -  (76) 

.T«  *  c^T 

R:  assuming  a  spherical  -fireball: 

R  -  (3  V  /  47f)  (77) 


Substituting  equation  (76)  for  V: 


R  - 


3  Y 


4Tr  <*'  +  c^T) 


1/3 


(78) 


A:  agaiOf  assuming  a  spherical  fireball: 


A  -  3  V  /  R 


(79) 


Substituting  equations  (76)  and  (78)  into  (79): 


A 


3  (477-/3)‘^^ 


Y 


+  c  T 

V 


2/3 


(BO) 
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Di-fferential  Eguation  -for  T (t) 

In  ordar  to  develop  a  first-order,  nonlinear,  ordinary 
differential  equation  for  T(t>,  it  is  noted  that  either  V  or 
T  can  be  used  as  the  independent  parameter.  Since  equation 
(32)  provides  dV/dt,  the  relations  developed  in  section  C 
above  can  be  used  to  convert  to  dT/dt. 

Applying  the  chain  rule  to  dV/dt: 


dV  _  dV  dT 
dt  “  dT  dt 


(81) 


which  can  be  solved  for  dT/dt.  First,  expressions  for  dV/dt 
and  dV/dT  in  terms  of  T  are  derived. 

Substituting  equations  (72),  (73),  (74),  and  (80)  into 

equation  (32): 


.  3_c 

dt  4 


_ Y _ 

aT^  ♦  c^T 


5/3 


(82) 


Differentiating  equation  (76)  with  respect  to  T: 


dV 

dT 


-  Y 


'^4aT^  ♦  c^T 


4  2 

(aT^  +  c^T)-' 


(83) 


Solving  equation  (81)  for  dT/dt  and  substituting  equations 
(82)  and  (83)  gives  the  desired  result: 


dT  _  -  3_c 

r.TT 

dt  4 

1 

1 

1 

1 

•<1 

1 

1 

1 

1 

1 

1 _ 

\^4aT^  +  c^jJ 
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E.  Initial  and  Final  Conditions 

Th»  initial  conditions  -for  the  one-temperature  model 
are  the  same  as  the  -final  conditions  for  the  three- 
temperature  model.  The  final  condition  for  the  one- 
temperature  model  is  hydroseparation,  which  occurs  when  the 
rate  of  growth  of  the  fireball  has  fallen  sufficiently  for 
the  acoustic  velocity  in  the  fireball  to  equal  dR/dt.  This 
occurs  at  about  300,000  C,  i.e.  about  25  eV  [Ref  1:653. 
Examples  computed  using  this  model  are  presented  in  the  next 
section.  The  endpoint  for  these  computations  is  taken  to  be 
T  ■>  25  eV.  More  detailed  consideration  of  the  conditions  for 
hydroseparation  is  a  question  of  hydrodynamics  and  is 
outside  the  scope  of  this  thesis. 
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VZ.  Nufflsrical  Results 


This  section  presents  graphs  of  the  temperatures » 
radii,  and  growth  rates  predicted  by  the  homogeneous  (three- 
temperature,  then  one-temperature)  model.  The  results  of  the 
early  time  model  are  shown  for  comparison.  The  homogeneous 
model  results  are  drawn  as  solid  lines  and 'the  early  time 
model  as  dashed  lines  in  figures  4  through  7,  below.  These 
results  are  for  a  20  kiloton  burst,  in  sea  level  air,  with  a 
square  yield  pulse  of  5.14493  shakes  duration.  The  curves 
labeled  "fit)"  are  included  to  indicate  the  shape  and 
duration  of  the  yield  pulse,  but  have  been  rescaled  in 
height  for  convenient  display,  and  thus  are  not  normalized 
to  unit  area  as  is  the  actual  f(t). 

Observations  regarding  the  important  features  of  these 
curves  are  made  in  this  section.  An  evaluation  of  the  model 
based  on  these  observations  is  presented  in  the  following 
section. 

A.  Material  and  Ihermal.  Radiation  Temperatures 

Figure  4  shows  T  and  T  ,,  for  20  shakes.  Important 

in  •tt 

features  of  this  graph  include: 

—  T^|. is  a  smoothly  varying  function  of  time 

—  T^  stays  at  T^  only  briefly,  but  is  held  out  of 
equilibrium  until  the  end  of  the  yield  pulse 
—  T^  >  at  all  times 

SB 


B.  R«diu« 

Figura  5  showm  th»  radius  as  a  ‘function  of  tins  for  20 
shakes  as  predicted  by  the  two  models.  Important  features 
include: 

—  The  homogeneous  model  radius  is  a  smoothly  varying 
function  of  time 

—  The  early  time  model  radius  converges  upon  the 
homogeneous  model  radius  within  about  S  shakes  after  the 
yield  pulse  ends 

—  The  models  are  very  similar  at  the  beginning  of  the 
yield  pulse. 

Growth  Rate 

Figure  6  shows  the  growth  rate,  dR/dt,  of  the  two 
models  for  20  shakes.  Important  features  of  the  homogeneous 
model  growth  rate  include: 

—  smoothly  varying  function  of  time 

—  monotone  decreasing 

—  never  exceeds  the  speed  of  light 

—  decreases  as  streaming  radiation  is  lost  after  the 
yield  pulse  ends. 

—  The  early  time  model  predicts  about  the  same  grotrth 
rate  once  the  spike  to  values  greater  than  c  decays. 
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Fig.  6  Growth  Rate  Comparison,  20KT 


Fig.  7  Hydrodynamic  Separation  Comparison  with  Nominal 
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D-  Hydrodynami c  Separation 

Hydrodynamic  separation  should  be  predicted  to  occur 
when  a  temperature  o-f  300,000^0,  i.e.  2S  eV,  is  reached,  at 

4 

a  time  of  about  0.1  millisecond,  i.e.  about  10  shakes,  for 
a  20KT  burst  CRef  1.651.  Figure  7  is  a  semi -log  graph  of 
for  both  models,  with  each  extended  until  it  predicts 
hydrodynamic  separaticsn.  The  nominal  case  described  above 
refers  to  a  total  yield  of  20  kilotons,  for  an  x-ray  yield 
of  about  15  kilotons.  The  curves  are  for  an  x-ray  yield  of 
20  kilotons,  but  the  difference  in  predicted  hydrodynamic 
separation  times  is  negligible  because  of  the  steep  decrease 
at  the  left  of  the  (semi-log)  curves. 

The  homogeneous  model  predicts  25  eV  at  15<X>  shakes, 
which  is  a  factor  of  6.7  underestimate  compared  to  10,000 
shakes.  The  early  time  model,  however,  predicts  364,000 
shakes,  which  is  a  factor  of  36.4  overestimate. 

Qj^bSC  Comparisons 

Appendix  A  presents  graphs  for  1  kiloton  and  4  megaton 
yields  and  for  various  pulse  shapes  and  durations,  and 
discusses  some  of  the  conclusions  which  may  be  drawn  from 
them.  The  features  seen  above  are  common  to  all  the  various 
yields  and  pulse  shapes. 


Th«  major  characteristics  of  the  model :Jbserved  from 
these  comparisons  are: 

— >  The  results  predicted  by  the  homogeneous  model  are 
very  insensitive  to  the  pulse  shape  and  duration*  so  long  as 
the  pulse  duration  is  short  compared  to  the  yield  (2  shakes 
at  IKT*  10  shakes  at  4nT) . 

—  The  results  predicted  by  the  early  time  model* 
specifically*  the  height  and  time  of  the  peaking  of  is 
sensitive  to  the  pulse  shape  and  duration*  even  for  short 
duration  pulses. 

Appendix  C  develops  a  modified  version  of  the 
homogeneous  model  which  accounts  for  spherical  divergence  of 
the  streaming  radiation.  Graphical  comparisons  to  both  the 
homogeneous  model  and  the  early  time  model  are  shown.  The 
significant  observation  drawn  from  those  comparisons  is  that 
the  modified  (divergent)  model  gives  essentially  the  same 
results  as  the  unmodified  (homogeneous)  model. 


VZI.  Evaluation  of  tha  Homogansoua  Modal 


This  aaction  avaluataa  tha  homogancKiua  modal  by 
conaidaration  o-f  thosa  -faaturea  which  wera  -found  (in  aaction 
XI)  to  ba  problama  with  tha  aarly  tima  modal.  Pradictiva  and 
padagogical  advantagaa  will  aach  ba  diacuaaad,  followad  by  a 
conaidaration  o-f  aoma  limitations  of  tha  modal*  and  o-f  its 
potantial  -for  growth  into  a  mora  sophisticated  modal. 

A.  ^adfctiye  Advantages  of  the  Homogeneous  Model 

1.  Ef -fective  Temperature.  *  smoothly  varying 

-function  of  tima*  as  would  reasonably  ba  expected.  Tha 
overall  shape  and  numerical  values  are  in  agreement  with  tha 
aarly  tima  model*  but  with  tha  peaking  (at  tha  transition 
from  burnout  to  diffusion)  smoothed  out.  This  is  a  result  of 
allowing  T  to  fall  below  T.  Independently  of  the  transition 

In  D 

from  streaming  to  diffusion*  and  of  accounting  for 
simultaneous  transport  by  streaming  and  diffusion  both. 

2.  Dacraasi ng  Fireball  Expansion  Rata.  dR/dt  is  a 
monotone  decreasing  function  of  tima  with*  in  particular*  a 
dacreasa  after  tha  yield  pulse  ends  due  to  loss  of  transport 
by  streaming  (which  is  four  times  mora  efficient  than 
transport  by  diffusion).  Tha  reasons  for  this  improvement 
are  tha  same  as  thosa  noted  in  tha  previous  paragraph. 


3.  Constraint  by  the  Speed  of  Light 


The  expansion 


rate  of  the  fireball  never  exceeds  the  speed  of  light, 
unlike  the  expansion  rate  predicted  by  the  early  time  model. 
The  reason  for  this  improvement  is  that  the  growth  rate  is 
modeled  by  consideration  of  the  conditions  at  the  boundary, 
using  a  local  energy  balance.  Since  energy  density  flows 
through  the  surface  at  a  maximum  speed  of  c,  and  the  energy 
density  accumulated  in  dV  is  greater  than  that  entering 
through  the  surface,  then  the  growth  rate  is  inevitably  less 
than  c. 

4.  T  versus  T  The  material  temperature,  T_, 
—m  — — — —  ~*eTT  m 

always  exceeds  the  effective  thermal  radiation  temperature, 
T^y  !■,  as  is  required  by  physical  considerations.  The 
homogeneous  model  is  successful  in  predicting  this  result 
because  it  does  not  decrease  E  by  dilution  (as  do  the 

iR 

moment  equations  of  the  early  time  model)  but  rather  it 
explicitly  models  the  conversion  of  E^  to  E^  by  absorption 
in  the  surface  absorption  layer  of  ambient  air. 

Pedagogical  Advantages 

!i!9  iQ^SSCSlS  5^fiC  Retarded  If  mg.  The  homogeneous 
model  avoids  the  confusion  of  retarded  time  and  the 
complexities  of  visualization  of  the  interaction  processes 
as  modeled  by  integrals  over  retarded  time. 
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2.  No  Moments  Equations.  The  homogeneous  model  avoids 
the  abstract  and  physically  unconvincing  use  o-f  moments 
equations  obtained  from  the  diffusion  equation. 

QSEiyid  by  Local,  and  Global  Energy  Balances.  The 
three-temperature  homogeneous  model  is  derived  using  only 
energy  balances  and  geometric  relationships.  Each  term  in 
the  balance  equations  can  be  interpreted  as'  describing  an 
individual  process  of  physical  significance.  This  approach 
engenders  the  sense  of  intuitive*  physical  understanding 
which  is  prerequisite  to  the  application  of  more 
sophisticated  techniques*  and  which  is  the  pedagogical 
objective. 

E*»i 1 V  Computerized.  The  homogeneous  model  leads 
directly  to  a  simple  finite  differencing  scheme 
(essentially*  a  forward  Euler  scheme).  This  scheme  can  be 
programmed  and  run  on  personal  computers  (which  is  becoming 
an  increasingly  important  resource  as  more  students  obtain 
them).  Doing  so  is  a  task  of  appropriate  size  to  be  assigned 
as  an  exercise  or  project. 

5.  Reasonable  N«.merical  Results.  The  results  obtained 


from  the  computer  exercise  suggested  above  are  numerically 
reasonable*  both  in  terms  of  curve  shapes  and  physical 
constraints*  and  in  terms  of  the  actual  numbers  computed. 


C.  Limitations  of  the  Homogeneous  Model 

The  accuracy  o-f  the  model  is  limited  by  that  of  the 
assumptions  made  in  its  development.  These  assumptions  are 
of  tiM)  kinds,  assumptions  about  the  physical  parameters  of 
the  air,  e.g.  heat  capacity  as  a  function  of  temperature, 
and  assumptions  about  the  distribution  and  transport  of 
energy  in  the  fireball. 

The  limitations  imposed  by  assumptions  about  the 
physical  parameters  are  not  of  great  conceptual 
significance.  If,  for  example,  the  assumption  of  constant 
heat  capacity  were  revised  to  an  assumption  of  some 
particular  temperature  dependence,  the  computcsr  program 
would  be  changed  only  by  the  introduction  of  a  subroutine  to 
evaluate  the  heat  capacity,  and  the  numerical  results  would 
change.  The  general  shape  of  the  resulting  curves  would  be 
unchanged,  as  long  as  a  reasonable  C^(T)  was  chosen,  and  the 
same  physical  constraints  (dR/dt  <  c,  etc.)  would  be  met. 

The  limitations  imposed  by  assumptions  about  the 
distribution  and  transport  of  energy  are  of  conceptual 
importance,  for  it  is  these  assumptions  which  constitute  the 
model.  The  assumption  of  uniform  spatial  distribution  of 
thermal  radiation  is  justified  by  the  fact  that  the  fireball 
is  optically  thin.  The  thermal  radiation  then  tends  to 
closely  couple  the  material  energy  in  different  regions  of 
the  fireball,  so  that  it  is  reasonable  to  assume  that  the 
material  temperature  is  uniformly  distributed. 
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Th»  assumption  that  ths  streaming  radiation  energy 
density,  E^,  is  also  uni-formly  distributed  in  the  -fireball 
must  be  recognized  as  the  major  conceptual  limitation  o-f  the 
■KJdel.  By  its  very  nature  as  streaming  radiation,  E^  must  be 
spherically  divergent,  and  so  cannot  be  uniformly 
distributed.  On  the  other  hand,  its  distribution  also 
depends  on  the  yield  rate  function,  f(t),  which  leads  to 
considerations  of  retarded  time.  The  reason  the  assumpticm 
was  made  was  to  eliminate  retarded  time  from  the  model.  The 
assumption  can  be  justified,  to  a  certain  extent,  by  the 
fact  that  it  results  in  reasonably  good  numerical 
predictions. 

The  assumption  of  spatially  uniform  E  can  actually  be 
replaced  by  taking  the  more  general  view  that  E  has  some 
(non-uniform)  distribution  within  the  fireball,  but  that  its 
interaction  rates  can  be  characterized  by  its  average  value 
within  the  fireball.  This  average  value  can  then  be  computed 
using  a  global  energy  balance.  The  important  limitation  is 
then  introduced  by  the  assumption  that  the  flux  of  streaming 
radiation  at  the  surface,  c  E^(R,t),  can  be  computed  by 
using  the  average  E^  as  an  approximation  for  the  surface 
value:  c  E  (R,t)  Z  c  E  (t)  From  this  viewpoint,  the 

S  9 

limitation  is  one  of  making  a  first-order  approximation 
rather  than  a  physically  unreasonable  assumption. 


68 


I 


Appendix  C  develops  a  modified  homogeneous  model  which 

assumes  that  the  streaming  radiation  is  spherically 

2 

divergent,  and  hence  distributed  as  1/r  .  This  ignores  the 
effects  of  absorption  and  emission  histories  on  the 
distribution.  The  resulting  approximation  is  that 

E.  <R,t)  :::  4Ej*''"^(t)  (ss) 

This  modification  causes  only  small  changes  in  the  predicted 
temperatures  and  radii,  indicating  that  the  homogeneous 
model  is  insensitive  to  the  distribution  of  streaming 
radiation,  whatever  it  may  be.  This  helps  explain  how  such  a 
simple  model  gives  good  numerical  results. 

Extension  to  More  Sgghisticated  Models 

Three  approaches  to  obtaining  better  numerical  results 
suggest  themselves:  better  input  data,  a  more  elegant  model, 
and  brute  force  application  of  the  model. 

Better  input  data  would  consist  of  using  detailed  data 
on  the  variation  of  the  physical  parameters  with 
temperature,  as  was  discussed  above. 

A  more  elegant  model  could  be  obtained  by  assuming  some 
non-uniform  spatial  distribution  for  one  or  more  of  the 
energy  densities,  and  using  it  to  obtain  a  higher  order 
approximation  for  the  surface  conditions  (which  determine 
dR/dt>  in  terms  of  the  average  conditions  (which  are 
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dBtermined  by  the  global  energy  balances).  Appendix  C 
presents  an  example  o-f  this  approach. 

The  homogeneous  model  presented  here  could  be  extended 
by  brute  force  to  include  radial  variations  in  temperatures 
and  energy  densities.  This  could  be  done  by  treating  the 
fireball  as  a  collection  of  concentric  spherical  shells. 
Each  shell  would  be  thin,  so  the  homogeneous  model  would 
apply.  The  boundaries  among  shells  could  be  treated  using 
the  same  methods  as  are  applied  at  the  surface  of  the 
fireball  in  the  homogeneous  model,  that  is  ,  the  same 
modeling  of  the  equilibration  and  transport  processes.  Such 
a  model  could  be  extended  past  hydrodynamic  separation  by 
treating  the  shells  as  Eulerian  coordinate  cells  for  hydro- 
dynamic  transport  simultaneously  with  the  radiative 
transport  model  developed  here. 
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VIII.  Conclusions  and  Recofflinandations 


A.  Pedagogical  Use  of  the  Homogeneous  Model 

The  goal  of  this  research  was  to  develop  and  evaluate  a 
model  of  the  radiative  growth  phase  of  fireball  growth, 
comparing  it  to  G.  C.  Pomraning's  early  time  model  CRef  21. 
The  model  sought  was  intended  for  pedagogical  use.  It  was 
therefore  to  be  based  on  simple  models  of  the  relevant 
physical  processes  and  to  provide  physically  and  numerically 
reasonable  results.  The  three-temperature  homogeneous  model, 
together  with  its  one-temperature  simplifications  (for 
application  after  thermal  equilibrium  is  closely  approached) 
has  been  developed  and  shown  to  be  such  a  model.  It  is 
recommended  for  pedagogical  use. 

B.  Computational  Use  of  the  Homogeneous  Model 

It  has  further  been  shown  that  this  homogeneous  model 
gives  numerical  results  which  are  more  reasonable  in  meetinr, 
necessary  physical  constraints  than  those  results  obtainf^d 
from  the  early  time  model.  Specifically,  the  growth  rate  of 
the  fireball  is  monotone  decreasing  and  always  less  than  thu 
speed  of  light,  as  predicted  by  the  homogeneous  model,  and 
the  effective  thermal  radiation  temperature  is  predicted  to 
be  a  smooth  function  of  time.  Overall,  however,  the 
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hf3iBog«n»ouB  model  predicts  temperatures  and  radii  Mhich  are 
quite  comparable  to  those  predicted  by  the  early  time  model. 
The  homogeneous  model  has  the  additional  advantage  of 
computational  simplicity  relative  to  the  early  tiaie  model. 
It  is  concluded  that  the  homogeneous  model  is  a  useful 
computational  method  of  at  least  first-order  accuracy.  It  is 
recommended  that  its  results  be  compared  with  those  of  a 
state-of-the-art  production  code  to  evaluate  its  actual 
accuracy  with  greater  confidence. 

C.  The  Homogeneous  Model  as  a  Basis  for  More  Sophisticated 
Models 

The  simple  models  of  radiative  emission,  absorption, 
and  transport  used  in  the  homogeneous  model  give 
surprisingly  good  numerical  results.  These  models  could  be 
applied  to  a  discrete  ordinates  method  to  obtain  more 
accuracy  in  the  radiative  stage  of  fireball  growth.  The 
discrete  ordinates  code  obtained  could  then  be  combined  with 
an  Eulerian  hydrodynamics  code  to  obtain  a  radiative 
hydrodynamics  model  for  the  later  stages  of  fireball  growth. 
These  approaches  are  recommended  for  consideration. 
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Appendix  A 

Comparison  of  Results  of  the  Two  Models 

This  appendix  presents  graphs  of  the  temperatures 
predicted  by  the  early  time  model  and  the  homogeneous  model. 
The  examples  presented  include  yields  of  one  kiloton  and 
four  megatons,  with  three  cases  of  yield  pulse  (shape  and 
duration)  for  each  yield. 

One  Kiloton  Examples 

Figure  A-1  shows  the  three  yield  pulses  considered  for 
the  one  kiloton  examples  and  indicates  the  choice  of  line  or 
dash  used  to  indicate  each  case  in  the  subsequent  figures. 
Figure  A-2  shows  the  effective  radiation  temperatures  for 
the  three  cases  as  predicted  by  the  early  time  model.  The 
temperatures  as  predicted  by  the  homogeneous  model  for  these 
same  cases  are  shown  in  figure  A-3  .  Conclusions  which  may 
be  drawn  from  these  figures  are: 

—  the  homogeneous  model  predicts  a  longer,  lower 
temperature  curve  for  a  longer,  lower  yield  pulse,  whereas 

—  the  early  time  model  predicts  a  longer,  higher 
temperature  curve  for  a  longer,  lower  yield  pulse) 

■—  the  homogeneous  model  predictions  are  insensitive 
to  yield  pulse  shape  as  long  as  the  pulse  is  short  (in 
comparison  to  the  total  yield),  whereas 

the  early  time  model  predicts  a  peaking  of  T^^^, 
the  timing  of  which  is  sensitive  to  pulse  shape. 


A-1 


Fig.  A-1  Yield  Pulses  for  1  Kiloton  EKseples 


Fig.  A-2 


Early  Time  Model  for  1  Kiloton  Yield 
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ft-3  Hofflogeneous  Model  for  1  Kiloton  Yield 
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Fig.  A-4  TemperetureB  for  Both  Models 
Mith  2  Shake,  Square  Pulse  of  1  Kiloton  Yield 
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T  iin©  <sh)a.kes> 


For  direct  comparimon  of  the  two  modele,  including  the 
material  temperature,  -figure  A-4  presents  both  T  and  T 

fli  I^T 

for  the  one  kiloton  yield.  The  yield  pulse  for  this  example 
is  case  2,  the  square  pulse  of  2  shakes  duration.  This 
comparison  demonstrates  the  good  general  agreement  of  the 
two  models. 


Four  Megaton  Examples 

Figure  A-5  shows  the  three  yield  pulses  considered  for 
the  four  megaton  examples  and  indicates  the  choice  of  line 
or  dash  used  to  indicate  each  case  in  the  subsequent 
figures.  Figure  A-6  shows  the  effective  radiation 
temperatures  for  the  three  cases  as  predicted  by  the  early 
time  modell  figure  A-7  is  an  expanded  view  of  this  same 
information.  The  temperatures  as  predicted  by  the 
homogeneous  model  for  these  same  cases  are  shown  in  figure 
A**8.  The  three  curves  in  this  figure  are  essentially 
identical.  Conclusions  which  may  be  drawn  from  these  figures 


—  the  homogeneous  model  is  very  insensitive  to  pulse 
shape  and  duration  for  this  large  yield,  whereas 

—  the  early  time  model  still  predicts  a  higher, 
longer  T^^ ^  pulse  for  a  longer,  lower  yield  pulse. 
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Fig.  A~5  Yi0ld  Pulses  for  4  Megeton  Examples 


Fig.  A-6  Early  Time  Model  T 
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for  4  Megaton  Yield 


Fig.  A-8  Homogsnvous  Mod»l  T  , .  for  4  nagaton  Yiald 


For  direct  comparison  o^  the  two  models,  including  the 


material  temperature,  figure  A-9  presents  both  T  and  T  ,, 

fli  St  t 

for  the  four  megaton  yield.  The  yield  pulse  for  this  example 
is  case  2,  the  square  pulse  of  5.14493  shakes  duration.  This 
comparison  again  demonstrates  the  good  general  agreement  of 
the  two  models. 


Fig.  A-9  Temperatures  for  Both  Models 
with  5.14493  Shake,  Square  Pulse  of  4  Megaton  Yield 
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Appendix  B 

Basic  Computer  Program 
•for 

Three-Temperature  Homogeneous  Model 

10  'Radiative  phase  •fireball  growth  program 

20  '  based  on  three-temperature  homogeneous  model 

30  *  by  LCDR  K.  A.  Mathews 

40  *  for  M.S.  Thesis  at  Air  Force  Institute  of  Technology 
50  WIDTH  LPRINT  80  'Set  Line  Printer  Line  Width 
60  DEFINT  I-N  'Set  Fortran-like  Integer  variable  naming 
70  DEF  FNCOEFKT)  « 

RHO.KAPPA.C  *  T.B/T.M)  *  (T.B/T.M)  *  (T.B/T.M) 
'function  call  for  temp,  variation  of  kappa 
80  DATA  8.9217E-2,  13.7,  400,  3,  3.14159,  4.96E10,  1.24EB 
90  READ  RHO.KAPPA.C,  A,  T.B,  C,  PI,  E.B,  C.V 

100  'Enter  problem  parameters  interactively 
110  HOMEIVTAB  2 

120  PRINT  "Enter  parameters  few  square  wave  yield  function:" 
130  PRINT 

140  HTAB  10:  INPUT  "Pulse  width  (shakes)  -  ",T. PULSE 
150  print:  PRINT 

160  PRINT  "Enter  Burst  Parameters:":  PRINT 

170  HTAB  10:  INPUT  "Yield  (kiloton)  -  ",  YIELD 

180  YIELD  *  YIELD  *  4.2E-t-12  'Convert  kilotons  to  Joules 

190  YDOT  »  YIELD  /  T. PULSE  'Constant  ydot  for  square  pulse 
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200  PRINT:  PRINT  "Enter  Quadrature  Parameters:": 


PRINT 


210  HTAB  lO:  INPUT  "Time  Increment  (shakes)  -  ",  DELTA. T» 
220  DT  *  DELTA. T#  indicates  double  precision  variable 

Time  is  carried  in  double  precision  -for  convenience  in 
assuring  output  times  round  to  desired  times  in 
printout 

230  print:  HTAB  10 

240  INPUT  "Number  of  steps  per  printout  -  ",  N. PRINT 
250  print:  HTAB  10 

260  INPUT  "Max  Time  to  which  to  compute:  T.quit  ",T.QUIT 
270  PRINT:  PRINT  "Commencing  computations" 

280  GOTO  400  ’to  initialize  variables 

300  ’Data  Lprint  Subroutine 
310  LPRINT  T,  V,  R,  R.DOT 

’Time,  Volume,  Radius,  dR/dt 
320  T.EFF  «  SQR  (SQR  (E.T  /  A)  ) 

’Computed  here  instead  of  in  main  routine  to  improve 
program  speed  of  execution;  T.eff  is  not  needed  for 
computations,  only  for  data  output 
330  LPRINT  E.S,  E.M,  E.T,  T.M,  T.EFF 
340  LPRINT 
350  RETURN 
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400  'Initialization  o-f  variables 
410  T  -  O:  T#  »  O 

420  E.S  ■  O  'will  be  reinitialized  by  special  first  time 

step,  as  will  E.fflt  all  are  set  to  zero  first 
for  benefit  of  initial  printout 

430  E.M  -  0 
440  E.T  *  O 
450  T.M  »  O 
460  T.EFF  -  O 
470  V  -  O 
490  R  -  O 
500  R.DOT  -  C 

600  'Initial  printout 

610  LPRINT  "Radiative  Phase  Fireball  Growth  Program:" 

620  LPRINT 

"  Version  4.3  —  Single  Precision,  No  Ground  Effects" 

630  LPRINT  "YIELD  <KT)  -  "» YIELD/4. 2E+12  'divide  to  get  KT 
640  LPRINT:  LPRINT  "Pulse  Shape  Parameters:" 

650  LPRINT 

"  Square  Wave  Pulse,  starts  at  t  >  0,  ends  at  t  «"|T. PULSE 
660  LPRINT  "  Yield  Rate  (J/shake)  -"lYDOT 
670  LPRINT:  LPRINT  "Quadrature  Parameters:" 

680  LPRINT  "  Time  Increment  (shake)  ■  "|DT 

690  LPRINT  "  Number  per  printout  -  "IN. PRINT 

700  LPRINT:  LPRINT  'space  printout  then  print  data  titles 


710  lprint:  lprint 


”T  (shaka)",*‘V  (n^3) ",  "Radius  (in>","dR/dt  (s/shaks)" 
720  LPRINT  "E.s  <J/«^3) ",  "E.a  <J/a'‘3)", 


•E.th  <J/a^3)","T.«  (aV)-,"T.»4f  (sV) " 


730  LPRINT 


740  GOSUB  300  ’Lprint  initial  conditions  undsr  ths  hsadings 


800  ’First  Time  Step:  inner  bdy  »  outer  bdy  at  t  «  O 
so  special  case.  Go  two  time  steps  at  once  using 
dR/dt  c  to  determine  volume,  and  volume  and  yield 
to  determine  E.s 

810  T  -  2  *  DT:  T#  -  2  *  DELTA. Ttt:  1-2 
820  T.M  -  T.B 

830  E.n  ■  E.B  ’Assume  Burnout  because  of  high  E.s  due  to 


small  volume 


840  R  -  C  *  T 


830  V-4*Pl*R*R*R/3 
860  IF  T<T. PULSE 

THEN  E.S  -  YDOT  *  T  /  V  -  E.M 
ELSE  E.S  -  YIELD  /  V  -  E.M 

’Yield  by  time  2*dt  is  ydot  *  t  unless  pulse  is 
already  over,  then  it  is  yield 
870  E.T  -  O  ’Assume  no  buildup  in  short  time  2«dt 
880  R.DOT  -  C 

900  GOSUB  300  ’Data  Printout  to  allow  manual  evaluation 

of  reasonableness  of  first  step,  then  start  main  loop 
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mm 


1000  'Main  Iteration  Routine 

1010  'Energy  Trans-format  ions  in  old  volume  V 

1020  IF  T.M  -  T.B 

THEN  RHO. KAPPA. C. OF. T  -  RHO. KAPPA. C 
ELSE  RHO. KAPPA. C. OF. T  »  FNC0EF(T.M> 

'Avoid  calculating  (T.b/T.b)  cubed  >  1  (-for  specKi) 

1030  DE.8M  •  RHO. KAPPA. C. OF. T  «  E.S  *  DT 

'This  is  potential  amount  to  be  absorbed,  it  will  be 
bounded  by  amount  required  to  maintain  burnout 
if  appropriate 

1040  DE.MT  =  RHO. KAPPA. C. OF. T  ♦ 

(A  *  T.M  *  T.M  *  T.M  *  T.M  -  E.T)  *  DT 
1030  IF  T.M  ->  T.B  AND  DE.MT  «<  DE. SM 
THEN  DE.SM  -  DE.MT 

'Bounding  if  at  burnout  to  prevef<t  raising  T.m  above  T.b 

1100  'Find  Volume  elements  for  radiative  transfer 

through  surface  in  dt.  The  energy  transferred  is 

that  present  in  a  shell  of  thickness  draK:*dt  (streaming) 

or  r«(c/4)«dt  (diffusion),  extending  inward  from  the 

surface  of  the  fireball.  If  the  shell  is  thick, 

its  volume  is  V  -  Vol (interior  to  shell) 

factually,  this  is  always  true>,  but  if  the  shell  is 

thin,  this  leads  to  errors  in  the  subtraction,  so 

its  volume  is  computed  as  4'ii’pi'»R'^2'»dr. 
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1110  IF  C  «  OT  <  .001  «  R  THEN  GOTO  1170 

1120  RS  -  R  -  C  *  DT 

1130  RT-R-C*DT/4 

1140  DOS  «  V  -  4  «  PI  *  RS  «  RS  «  RS  /  3 

1130  DVT  -  V  -  4  ♦  PI  *  RT  *  RT  *  RT  /  3 

1160  GOTO  1200 

1170  DVS  -  4*P1*R*R*C*DT 
1180  DVT  ■PI*R*R*C*DT 

1200  'Assume  streaming  only  causes  dR/dt, 
compare  to  c/4  to  verify  assumption 
1210  DV  -  DVS  «  E.S  /  (E.S  *  E.M) 

1220  V.NXT  -  V  +  DV 

1230  R.NXT  -  (.73  ♦  V.NXT  /  PI)  (1/3) 

1240  R.DOT  ■  (R.NXT  -  R)  /  DT 
1230  IF  R.DOT  >-  C/4 

THEN  DVT  -  DV:  GOTO  1400 
'Assumption  correct,  set  dvt  to  dv  so  thermal 
does  not  contribute  to  E.m  in  growth  volume, 
but  E.t  is  decreased  by  dilution 


1300  'Otherwise,  recompute  v.nxt,  etc.,  to  include  thermal 
diffusion  in  growth,  and  recompute  needed  parameters 
1310  DV  ■  (DVS  ♦  E.S  +  DVT  *  E.T)  /  (E.8  +  E.T  +  E.H) 

1320  V.NXT  «  V  +  DV 

1330  R.NXT  -  <.75  *  V.NXT  /  PI)  Cl/3) 

1340  R.DOT  -  (R.NXT  -  R)  /  DT 

1400  'Energy  balances  on  new  volume  for  each* of  the  energies 
1420  DILUTED  >  (V-^DV/2)  /  V.NXT 
1430  IF  T>T. PULSE  THEN  DE. 8. YIELD  -  O 
ELSE  IF  T+DT  <  T. PULSE 

THEN  DE. 8. YIELD  -  YDOT  *  DT  /  V.NXT 
ELSE  DE.S. YIELD  -  YDOT  *  CT. PULSE  -  T)  /  V.NXT 
'Compute  addition  to  E.s  based  on  yield  delivered  during  dt 
1440  E.S  •  E.S  *  (1  -  DVS/V.NXT) 

-  DE.Sn  «  DILUTED  DE.S. YIELD 

1450  E.M  -  E.H  (DE.Sn  -  DE.NT)  *  DILUTED 
1460  E.T  -  E.T  ♦  <1  -  DVT/V.NXT)  +  DE.HT  *  DILUTED 

1500  'Set  variables  to  next  values 

1520  V  -  V.NXT 

1530  R  -  R.NXT 

1540  T#  -  T#  +  DELTA. T# 

1550  T  -  T# 

1560  T.H  -  E.M  /  C.V 
1570  I  -  I  ♦  1 
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1600  'End-o-f-PhasB  Criteria  checks 


1610  IF  T.M  <■*  25  OR  T  >■  T.QUIT 

THEN  60SUB  300:  GOTO  1700  'Line  Print  final  data 
1620  ’Printout  Criteria  Check 
1630  IF  (I  MOD  N. PRINT)  -  O 
THEN  I  -  O:  GOSUB  300 

’reset  I  to  zero  to  avoid  overflow  of  integer  data  type 
1640  GOTO  1000  ’Loop  to  next  iteration 

1700  ’End  of  Radiative  only  fireball  growth. 

Hydroscparation  assuffled  to  occur  at  temp  of  25  eV 
unless  T.quit  is  reached  first 
1710  LPRINT  CHR4(12)  ’Form  Feed 

1720  GOTO  110  ’Loop  to  parameter  input  routine  for  next 
problem. 

Notes: 

Everything  on  a  numbered  line  subsequent  to  a  ’  is 
a  remark  and  is  ignored  by  the  interpreter  and  compiler. 

—  Multiple  executable  statements  on  a  single  numbered 
line  are  separated  by  colons 

—  LPRINT  prints  on  the  line  printer!  PRINT  writes  to 
the  display  screen 

<—  This  listing  follows  the  syntax  of  BASIC-GOp 
revision  5.2  for  the  Apple  3C  Computer  with  Z-GO  auxiliary 
processor p  and  is  compatible  with  both  the  interpreter  and 
coflipiler  for  that  language,  as  distributed  by  Microsoft. 
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—  stability  of  the  Algorithm; 

This  algorithm  is  only  conditionally  stable. 
As  equilibrium  is  approached,  the  time  increment,  dt,  must 
be  decreased  to  maintain  stability.  However,  as  dt  is 
decreased,  the  subtractions  in  lines  1440  through  1460 
introduce  rounding  error  inaccuracies.  Increasing  the 
floating  point  word  length,  such  as  by  using  double 
precision,  only  delays  the  loss  of  accuracy,  allowing  a 
shorter  dt.  The  solution  to  the  difficulty  is  to  use  this 
program  to  the  point  of  instability,  and  inspect  the  results 
to  find  reliable  data  on  the  approach  to  equilibrium.  These 
results  are  then  used  as  the  initial  conditions  for  the 
numerical  integration  of  the  one-temperature  differential 
equation  for  T(t>.  Use  of  single  precision  (32  bit)  floating 
point  arithmetic  for  the  combined  model  readily  provides 


numerical 

results 

with 

enough 

significant 

digits  for  a 

tolerance 

of  1 ess 

than 

1  eV. 

This  method 

was  used  in 

computing  the  examples  presented  in  this  thesis. 
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Appendix  C 

'  Modification  of  the  Homogeneous  Model 

to  Include 

Spherical  Divergence  gf  the  Streaming  Radiation 

Objectives 

The  objective  of  the  main  body  of  this  thesis  has  been 

to  present  the  simplest  possible  model  for  the  radiative 

growth  phase  of  fireball  formation  in  an  atmospheric  nuclear 

burst.  As  discussed  in  section  VI I. C,  a  significant 

limitation  of  the  homogeneous  model  is  its  assumption  of  a 

uniform  density  of  streaming  radiation.  This  assumption 

2 

seems  unreasonable  because  of  the  1/r  radial  attenuation 

caused  by  the  spherical  divergence  implicit  in  the  nature  of 

the  streaming  radiation  (as  streaming  from  a  point  source). 

This  appendix  will  develop  the  modifications  to  the 

2 

homogeneous  model  to  include  this  1/r  radial  attenuation. 

Factors  Affecting  Radial  Dlstri.buti.on  of 

There  are  actually  three  factors  which  affect  the 
radial  distribution  of  streaming  radiation  at  any  fixed 
instant  of  tiffle,  t,  of  which  divergence  is  only  one.  The 
other  two  are  the  yield  history  which  resulted  in  the 
emission  of  streaming  radiation  present  at  time  t  and  the 
absorption  to  which  the  radiation  at  each  radius  has  been 
subjected  in  reaching  that  radius  by  time  t.  To  include 


J 
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thM»  lattsr  two  af facts  would  rsquirs  axplicit  uss  of 


rstardsd  tiins  and  integrations  over  retarded  time,  as  is 
done  in  the  early  time  model.  The  objective  of  this 
appendix,  however,  will  be  to  isolate  the  effects  of 
divergence  from  those  of  emission  history  and  absorption 
history,  and  hence  avoid  the  use  of  retarded  time. 

Assumptions  of  the  Divergent  Model 

In  avoiding  the  use  of  retarded  time,  the  divergent 
model  (as  the  modified  model  will  be  called)  assumes  that 
absorption  is  negligible  with  respect  to  divergence  in 
determining  the  radial  distribution  of  the  streaming 
radiation,  and  that  the  yield  rate  has  been  approximately 
constant  over  the  interval  in  which  the  streaming  radiation 
present  at  any  time  t  was  emitted.  These  conditions  are 
effectively  met  during  the  emission  period  of  a  square- 
shaped  yield  pulse  of  short  duration  relative  to  the  total 
yield.  For  this  reason,  the  divergent  model  will  be 
evaluated  for  three  example  cases,  a  one  kiloton  burst  with 
a  2  shake  square  pulse,  and  twenty  kiloton  and  four  megaton 
bursts  with  a  S. 14493  shake  square  pulse.  These  choices  also 
allow  comparison  with  the  examples  presented  in  sections  ZI 
and  VI  and  in  Appendix  A. 
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Method 


The  vieNpoint  suggested  in  section  VII. C  will  be  taken. 
Thus,  the  interactions  o-f  material  and  radiations  inside  the 
fireball  are  assumed  to  be  characterized  by  the  average 
values  of  the  energy  densities,  E_,  and  £«..  The 

S  PI  b 

interactions  at  the  surface  (in  the  thin  absorption  layer  at 

the  surface)  of  the  fireball  will  be  determined  by  the 

average  values  for  the  material  and  thermal  mergies,  but  by 

the  surface  value  for  streaming  radiation,  E^(R).  This 

surface  value  is  proportional  to  the  average  value  with  a 

proportionality  constant  which  is  to  be  determined  by  the 

2 

assumption  that  E^<r)  ~ 1/r  . 


Derivation  of  E^<R)  as  a  Function  of  E^ 

2 

Using  the  assumption  of  1/r  attenuation,  let  the 
constant  B  (possibly  a  function  of  time)  be  defined  by: 

E^(r)  »  B  /  r^  (C-1) 

Then,  the  average  value  of  the  streaming  radiation  density 
is  evaluated  by: 


(C-2) 


TABLE  C-1 


Divergent  Three-Temperature  Model  Balance  Eguatione 


TABLE  C~2 


SlySCflSDt  IbCg®  Temperature  Model  Differential  Equationa 


A  -  4  77  R" 


A  R  /  3 


5t 


(VE  ) 
dt  m 


-  V  dE  /dt  +  <E  /3) 

am  a  dt 


+  Y  f(t)  -  c  (E  /3)  A 


V  CdE  /dt  -  dE  ./dtl 
am  fflt 


E 

m  dt 


<C-14) 


(C-15) 


(C-16> 


If  dR/dt  >  c/4  then 


ZZ  <VE.  )  «  V  C<T  )  CaT^  -  E*.3 

Ou  u  in  fn  b 


and 


dV 

dt 


A  c  <E^/3)  /  C<E^/3)  +  E„3 
a  am 


(C-17) 


(C-18) 


OtherMi ae 


5t  -  ''  at  *  I  ^  « 


and 


(C-19> 


“  ■  A  c  C<E  /3)  +  E./43  /  C  (E  /3>  +  E.  +  E  3 

db  a  c  a  t:  n 


Where: 


and 


dE_./dt  -  C<T  )  CaT.-E.l 
mt  m  me 


If  T  »T.  and  CaT  -E.  3  <  E 
mo  me  a 


(C-20) 


(C-22) 


Then 


dE  /dt  -  C(T^)  CaT>E.3 

sm  m  fli  b 


<C-23) 


Otherwiae  dE  /dt  ■  C(T  )  E 
am  n  i 
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Computer  Program  Modifications 

The  computer  program  listed  i«  appendix  A  was  modified 
to  include  the  effects  of  spherical  divergence  of  the 
streaming  radiation.  The  modifications  to  the  program  follow 
directly  from  the  modifications  to  the  balance  equations 
described  above.  The  altered  program  lines  are  listed  below: 

1210  DV  >  DVS  *  (E.S/3)  /  (E.S/3  E.M) 

1310  DV  =  (DVS  *  (E.S/3)  +  DVT  ♦  E.T) 

/  (E.S/3  +  E.T  +  E.M) 

1440  E.S  «  E.S  «  (1  -  ( (DVS+2*DV)/3)/V.NXT) 

-  DE.SM  «  DILUTED  DE.S. YIELD 

The  global  energy  balance  (upon  which  Line  1440  is  based) 
has  been  rearranged  to  improve  numerical  performance. 

Results  of  the  Divergent  Modef 

Figure  C-1  shows  the  material  and  thermal  radiation 
temperatures  predicted  by  the  divergent  model  for  a  one 

kiloton  yield  with  a  2  shake  square  pulse.  Figure  C->2 

compares  T^^  ^  from  figure  C-1  with  T^^^  as  predicted  by  the 
homogeneous  model  and  the  early  time  model  for  the  same 

yield  pulse.  Figures  C-3  and  C-4  present  this  same 

information  for  the  case  of  a  four  megaton  yield  with  a 
5.14493  shake  square  pulse. 


Fig.  C-1  Oivargent  Mod»l  Tsmperatures 
for  One  Kiloton  Yield  Nith  2  Shake  Square  Pulse 


Fig.  C-2  Effective  Temperatures  for  the  Three  Models 
for  One  Kiloton  Yield  with  2  Shake  Square  Pulse 
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Fig.  C-3  Divergent  Model  Temperatures 
■for  Four  Megaton  Yield  with  5.14493  Shake  Square  Pulse 


Fig.  C-4  E-ffective  Temperatures  for  the  Three  Models 
for  Four  Megaton  Yield  with  5.14493  Shake  Square  Puls 


Figure  C-S  compares  the  radius  predicted  by  the 
divergent  model  with  that  predicted  by  the  early  time  model 
and  by  the  (unmodified)  homogeneous  model  for  the  case  of 
twenty  kiloton  yield  with  a  5.14493  shake  square  pulse.  This 
is  the  case  used  for  the  examples  in  sections  II  and  VI. 


Fig.  C-5.  Radius  for  the  Three  Models 
for  Twenty  Kiloton  Yield  with  a  5.14493  Shake  Square  Pulse 

Observations 

The  temperature  curves  predicted  by  the  divergent 
model  do  not  differ  qualitatively  from  those  of  the 
(unmodified)  homogeneous  model.  The  addition  of  spherical 
divergence  does  result  in  a  somewhat  improved  agreement  with 
the  early  time  models  but  the  use  of  a  state-of-the-art 
computer  code  to  determine  the  "right"  answer  would  be 
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nMdsd  to  choose  between  the  two  versions  of  the  honogeneous 
model • 

At  -first  thought,  it  might  seem  that  reducing  the 
streaming  radiation  by  a  factor  of  three  would  decrease  the 
growth  rate  by  a  factor  of  about  three,  but  figure  C-S  shows 
this  not  to  be  the  case.  There  are  several  reasons  for  this 
insensitivity: 

—  If  E  is  quite  large  compared  to  E^,  then  even  E/Z 
is  large  compared  to  E^,  and  dR/dt  is  near  the  speed  of 
light  and  is  insensitive  to  E^. 

*-*-  If  E^  is  large  compared  to  E^,  then  E^  determines 
dR/dt,  rather  than  E^. 

—  If  dR/dt  is  decreased  somewhat  (which  it  is),  then 
the  volume  wii:  ^'e  smaller  (by  the  ratio  of  radii  cubed)  and 
the  radiation  energy  will  have  a  higher  average  density, 
which  tends  to  offset  the  effect  of  the  1/3  factor  in  E^(R), 
and  prevent  dR/dt  from  decreasing  as  much  as  anticipated. 

In  any  event,  the  divergent  model  provides  very  good 
agreement  (in  predicting  the  radius)  with  the  early  time 
model  during  the  early  part  of  the  yield  pulse  for  a  square 
pulse.  This  is  as  sKpected,  for  during  this  period, 
absorption  is  small  compared  with  divergence,  and  the 
streaming  radiation  yield  history  is  constant. 
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Conclusions 


1.  Tho  inclusion  o-f  spherical  divergence  of  the 
streaning  radiation  provides  some  improvement  in  agreement 
with  the  early  time  model  with  respect  to  temperatures,  and 
more  improvement  with  respect  to  radius. 

2.  The  improvement  is  quantitative  and  small,  rather 
than  qualitative.  It  is  doubtful  whether  the  additional 
complications  are  justified  by  the  imprbvement  if  the 
homogeneous  model  is  used  pedagogical ly.  Some  students  may 
object  to  the  assumption  of  uniform  streaming  energy  density 
in  the  homogeneous  model.  The  modifications  developed  in 
this  appendix  could  be  of  value  as  a  response  to  such 
objections. 

3.  The  homogeneous  model  has  been  shown  to  be  very 
insensitive  to  large  variations  in  the  assumed  distribution 
of  the  streaming  radiation.  This  tends  to  explain  the 
reasonably  good  numerical  results  obtained  by  the  model,  and 
is  the  most  important  conclusion  of  this  appendix.  The 
assumption  of  spatially  uniform  may  be  a  poor  assumption 
in  and  of  itself,  but  the  model  is  so  insensitive  to  the 
distribution,  whatever  it  may  be,  that  the  assumption  is 
adequate  for  the  purpose. 
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Appendix  D 

M»thod  for  Recurslye  ^valuation  of  the  Integrals 
in  the  Earl.y  Hme  Model  lntgQ^o~Pi'f'fBrential  Eouatione 

Stateeent  of  the  Problem 

Pomraning's  early  time  model  CRef  23  consists  of  a 
linked  system  of  integro-di'f'ferential  equations,  for  the 
burnout  phase.  In  reference  3,  Pomraning  presents  a  finite 
differencing  scheme  for  integrating  this  system.  This  scheme 
requires  a  trapezoidal  integration  of  each  of  two  integrals 
(which  are  coefficients  in  the  differential  equations)  at 
each  time  step.  This  approach,  although  straight-forward, 
requires  maintaining  arrays  of  values  of  the  variables  in 
the  integrands,  amounting  to  thousands  of  floating  point 
numbers,  as  well  as  causing  expense  in  execution  time.  These 
drawbacks  would  have  precluded  implementation  of  the  model 
on  the  personal  microcomputer  used  for  this  research.  It  was 
therefore  necessary  to  derive  a  recursive  system  for 
evaluation  of  these  integrals.  This  appendix  presents  the 
recursive  method  which  was  used. 

Integrals  to  be  Evaluated 

The  integrals  to  be  evaluated  appear  in  equations  (41) 
and  (42)  of  reference  2  (equations  (42)  and  (43)  of 
reference  3).  Using  T as  the  dummy  of  integration,  rather 
than  Pomraning's  t',  the  integrals  to  be  evaluated  are: 


i<t)  -  r  cv<r)/v(t)3  dr  <D-i> 

'0 


and 


H(t) 


/J 


Sit  -  -“>  dr  /  (D-2) 


Mhara 


-  yOC 


(D-3) 


and 


S<t)  -  b  CaT^’  -  E^(t)3 


(D-4) 


Racursion  Relation  -for  IIH 

Let  indicate  I(t^),  where  t^  is  the  i'th  time  in  a 
■finite  differencing  scheme,  and  assume  that  (and  the 

other  variables  at  time  t^)  are  known,  as  is 
objective  is  to  find  a  formula  for  in  terms  of  and 

Let  the  time  increment  be  jt,  then 


“  ^i+1  ”  ^i 


(D-5) 


Using  this  notation,  equation  (D-l>  evaluated  at  t-i-it 
becomes: 


<D-6> 


Segmenting  the  interval  of  integration  and  rearranging: 
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1 


*i+l  “  ^-b<t-T)  ^-b-gt 

♦  ^^CV<T)/V^^j3  ^-b<t-r+gt) 

The  first  integral  in  equation  <D-7)  is  proportional  to 
and  the  second  can  be  evaluated  using  the  trapezoidal  rule: 


[  »  (V  /V  )  1 

^i+l  ^i  i+1  i 

+  Cl  +  e"**'^*!  gt  /  2 


(D-8) 


For  cofflputaftonal  efficiency,  using  a  constant  tins 

increment  gt,  let 

TJ  -  e"**'^^  (D-9> 

and 

Cj  -  yo  c  K  a  (D~10) 

Then  E^<t)  *  I(t),  and  Pomrtming's  equation  (41)  becomes 

the  desired  recursion  relation: 


-  T}^i  ♦  t  <l  *7;>  /  2  CD-11) 

where  indicates  E^(t^). 

This  is  the  basis  for  line  lOSO  of  the  program  listed 
in  Appendix  E. 
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Conversion  to  Retarded  Ti.flje  Integration  for  Hit) 

It  is  convenient  to  maintain  a  list,  of  S(t^),  but 
awkward  to  convert  these  values  to  correspond  to  steps  of  dr 
for  S(t-<R-r>/c) .  Thus,  it  is  efficient  to  change  variables 
-from  r  to  r,  where  T  is  defined  by 


where 


c  (r  -  t') 


(D-12) 


t'  -  t  -  R(t)/c 


(D-13) 


The  physical  significance  of  t’  is  that  t*  is  the  time  of 
emission  of  the  streaming  radiation  which  reaches  the 
surface  of  the  fireball,  R(t),  at  time  t.  Similarly,  T  is 
the  time  of  emission  of  the  streaming  radiation  which 
reaches  radius  r  at  time  t.  Thus,  these  are  retarded  times 
in  the  sense  of  the  retarded  potentials  of  electromagnetic 
theory. 

The  transformed  integral  for  H(t>  is: 

r  ^  2 

H<t)  -  /  (r-  t')'^  S(r>  dr  (D-14) 

't' 


Approach  for  H^t ). 

The  related  integral  forms,  F(t>  and  G(t>,  will  be 
needed  for  the  recursion  relation  for  H(t>.  These  are 
defined  by: 
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r  t 

F(t)  - 

/  Sir)  dr 
» 

iD-19) 

b 

Bit)  - 

1  ir-  t*)  Sir)  dr 

iD-16) 

It  nay  bn  notnd  that  these  forms  are  similar  to  Vol terra 
integrals  with  difference  kernels.  However,  t*  is  an 
implicit  function  of  t,  which  complicates  use  of  Leibnitz* 
rule  such  that  no  useful  simplifications  could  be  found  in 
that  direction. 

Techniques  similar  to  those  used  for  l(t>,  above,  lead 
to  a  recursion  relation  for  H^^^  in  terms  of  H^,  6^,  and 
other  known  quantities.  A  recursion  relation  is  then  needed 
for  And  is  found  to  be  in  terms  of  6^,  ,  etc.  The 
recursion  relation  then  needed  for  i*  only  in  terms  of 
F^  and  known  quantities,  completing  the  recursive  system.  To 
simplify  the  presentation,  the  recursion  relations  are 
derived  in  increasing  order  of  complexity,  from  F  to  H. 

Recur si on  Relation  for  F (t) 

Segmenting  the  interval  of  integration  of  equation  (D- 
19),  and  evaluating  at  for 


S<r)  dr 


S<r)  dr 


<D-17) 
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Tha  intagrals  in  equation  (D-17>  are  evaluated  using  the 
trapezoidal  rule,  which  produces  the  desired  recursion 
rel ation: 


whore 


i+1 


St’  (sr^.  *  s,’>  /  2 

1+1  i 


(D-18) 


S'  -  S(tJ>  (D-19) 

Sf  -  tJ^j  -  t'  -  St  -  fiR/c  (D-20) 

5R  -  -  Rj  <D-21) 


and  where  *  ^^^i+1^  obtained  by  linear 
interpolation  in  the  list  of  S^.  (Note  that  so 
it  is  always  possible  to  perfore  the  interpolation.)  8*^^  is 
retained  for  use  as  8'  at  the  next  tiee  step,  so  only  one 
(interpolated)  table  look-up  is  required  at  each  time  step. 


Recursion  Relaticxi  for  B(t) 

Development  of  the  recursion  relation  for  will  be 

simplified  by  the  following  notation:  define  6(a,b)  by 


8(a,b)  -  f  (r-  t'  -^t')  S(r)  dr  (D-22) 

'a 

Then,  8^^^  -  8(tJ^j,tJ)  -  8(t'+St» , t+fit)  where  the 

subscript  "i*  is  dropped.  Using  this  notation  and  segmenting 

0-6 


th«  Interval  o-f  integration  as  be-fore 


-  G<t',t)  +  B(t,t+Jt)  -  B(t»,t'+{t')  <D-23> 


The  three  terms  on  the  right  o-f  this  equation  will  now  be 
evaluated.  Note  that  G(t',t)  ^  because  o-f  the  -ft'  in  the 
integrand  o-f  equation  <D-22).  Since  St'  is  a  constant  with 
respect  to  the  integration,  however,  the  integration  can  be 
simpli-fied  as  -follows: 


G(t',t) 


<^-  t')  S(r)  dr 


-  ft 


’C 


S(r)  dr 


(D-24> 


Recalling  equations  <0-lS)  and  (D-16) : 


B<t',t)  «  Gj  -  ft' 


<D-25) 


The  second  term  in  <D-23)  is  evaluated  by  the 


trapezoidal  rule  as: 

G(t,t-*^ft)  -  ft  C<t+ft-t'-ft') 

+  (t-t'-St'>  S^3  /  2 


(D-26) 


This  is  simplified  by  substituting  the  definition  of  the 
retarded  time,  ■  ^i+1  ”  that 

t+ft-t'-St'  -  tj^j  -  tJ^j  -  (D-27) 
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ii«hich  leads  to 

t-t'-St'  »  -  5t  (D-28> 

Substituting  equations  (0-27)  and  (0-28)  into  (0-26)  gives 

B(t,t+5t)  -  6t  ®i  ‘“^i+i^^  -«rt)3  /  2  (0-29) 

The  third  term  in  equation  (0-23)  is  evaluated 
sifflilarly  as: 


G(t',t'+5t')  -  5t»  c(t'+it'-t'-^t')  sr  . 

'  ‘  (0-30) 

+  (t'-t'-it')  Sri  /  2 
1 

which  simplifies  algebraically  to 

G(t',t'+5t’‘ )  »  -  8t'^  S'  /  2  (0-31) 

Substituting  equations  (0-25),  (0-29)  and  (0-31)  into 

(0-23)  gives  the  desired  recursion  relation  for  6(t): 


G 


i+1 


sr  /  2 

1  ''  1  1 
+  St  CSj^^  ♦  Si  -  5t)3  /  2 


(0-32) 


This  recursion  relation  is  entirely  in  terms  of  known 
quantities  and  requires  no  additional  table  look-ups. 
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Rscursion  Relation  -for  H<t> 


The  approach  -for  H<t)  is  the  same  as  that  -for  Git).  Let 
Hiayb)  be  de-fined  by: 

r**  2 

H<a,b)  «  j  <r-  t'  -  Sir)  dr  (D-33> 

‘'a 

Then  equation  (0-14)  evaluated  at  t^^^  becomes 

”l+l  "  H(t,t+5t)  -  H<t»,t»+^t')  (D-34) 

The  integral  -for  H(t*«t)  is  evaluated  by  expanding  the 
square  in  the  integrand  and  expressing  the  result  as  the  sum 
o-f  three  integrals,  as  -follows: 

H<t',t)  -  f  ir-t')^  Sir)  dr 

''t' 

-  2  ft'  f  ir-  t*)  Sir)  dr 
2  r* 

+  /  Sir)  dr  iD-35) 

-'t' 

These  expressions  are  the  reason  for  the  original  choice  of 
the  definitions  of  Fit)  and  Git).  Substitution  of  those 
definitions  gives 

Hit»,t)  -  -  2  ^t'  Gj  +  ^t'^  Fj  iD-36) 
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Th»  macond  term  in  equation  (D-34)  ie  evaluated  using 


the  tr«pszoidal  rule  to  obtain: 


H(t,t-^ft>  -  St  C<t+St-t'-«t')^ 

+  (t-t'-St')^  Sj3  /  2 


(D-37) 


Substituting  equations  <D-27>  and  (D-2B)  gives  the 
computationally  convenient  form: 


H<t,t+St)  -  St  C<R^^j/c)^ 

+  -  St)^  Sj3  /  2 


(0-38) 


The  third  term  in  equation  (0-34)  is  evaluated  using 
the  trapezoidal  rule  to  obtain: 


H<t',t’+StM 


St'  C(t'+St'-t'-St')^ 

+  Ct'-t'-St')^  SJ3  /  2 


(0-39) 


which  simplifies  algebraically  to: 

H(t',t'+St')  -  St'^  SJ  /  2  (0-40) 

Substituting  equations  (0-36),  (0-38)  and  (0-40)  in  equation 
(0-34)  produces  the  desired  recursion  relation  for  H(t): 


H 


i+1 


H. 


-  2  St'  G, 


St 


,2 


-  St 


,3 


S'  /  2 


+  St  C(Rj^j/c)^  Sj^j  +  Sjl  /  2 


(0-41) 
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TABLE  D-1 


Recurs! ve  System  for  Evaluation  of 
jthe  Retarded  lime  Integrmli  H<t) 


^i+1 

-  Fj  +  5t  +  Sj)  /  2 

-  ^t'  <S'^j  +  SJ)  /  2 

(D-18) 

®i+l 

»  -  St*  Fj  +  ^t'^  S'  /  2 

+  St  R.^j/c  +  S.  (R.^j/c  -St)2  /  2 

CD-32) 

»i+l 

-  Hj  -  2  St*  6.  +  5t'^  Fj  -  ^t'^  SJ  /  2 

+  St  C<Rj^j/c)^  ~  SjJ  /  2 

(D-41> 

Nhere: 

t'  -  t  -  R(t>/c 

(D-13> 

SJ  »  S(t') 

(D-19) 

^t'  -  t'^j  -  t'  -  Jt  -  SR/c 

(D-20) 

SR  -  Rj^j  -  Rj 

(D-21) 

and 

where  SJ_^j  -  S(tJ^j)  is  obtained  by 

linear 

interpolation  in  the  list  of  S^. 

Conclusions 


Ths  recursive  system  -for  evaluation  o-f  HCt)  is 
summarized  in  Table  0-1.  The  system  has  the  advantages  of 
being  entirely  algebraic  and  of  requiring  only  one 
(interpolated)  table  look— up  at  each  time  step.  If  a 
constant  time  incrementf  jt,  is  usedp  the  table  look-up 
becomes  an  algebraic  calculation,  rather  than  a  search. 
These  features  speed  program  execution  significantly  (from 
80  hours  using  the  trapezoidal  integrations,  estimated  based 
on  the  execution  time  of  the  example  case  described  below, 
to  12  minutes  using  the  recursion  relations).  They  also 
minimize  computer  storage  space  requirements,  since  only  a 
list  of  S  values  need  be  maintained.  This  recursive  system 
is  the  basis  of  lines  1140  through  1230  of  the  program 
listed  in  Appendix  E. 

Benchmark  Comparison 

The  numerical  accuracy  of  this  method  was  verified  by 
comparison  with  computation  directly  by  trapezoidal  rule. 
The  example  case  used  was  evaluation  at  a  single  time  step 
at  the  end  of  the  burnout  phase  for  a  four  megaton  burst. 
The  results  of  the  two  calculations  agreed  identically  (to 
the  limit  of  precision  carried).  To  verify  the  entire 
implementation  of  Pomraning's  model,  a  more  complex  program, 
of  which  the  program  listed  in  Appendix  E  is  a  subset,  was 
used.  This  larger  program  included  the  effects  of  the 
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ground.  Tho  oxamplo  choson  (4  MT  burst  at  1  m  hsight, 


exponential 

pul se  shape 

with  ce  -  ^  - 

0.6 

shake 

duplicated 

that  used  i n 

reference  2.  This 

example  case 

served  as 

a  benchmark. 

The  results,  T  (t) 

in 

and 

T.ff<t). 

obtained  with  this  program  matched  those  graphed  in 
reference  2,  within  the  precision  with  which  that  graph 
could  be  read.  This  agreement  was  obtained  -for  the  burnout 
and  diffusion  phases,  verifying  the  validity  of  both 
programs. 
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Appendix  E 

Be«ic  Computer  Program 
•for 

Burnout  Phase  ot  the  Earl,j^  liES  Model 


This  version  of  the  program  assumes  a  triangular  yield 
pulse  of  isosceles  shape.  It  is  a  more-'or-less 
straightforward  implementation  of  Pomraning's  finite 
differencing  scheme  CRef  33.  fWi  important  exception  is  that 
the  explicit  computation  of  integrals  as  coefficients  at 
each  time  step  has  been  replaced  by  recursive  schemesy  as 
derived  in  appendix  D.  Variable  names  follow  the  notation  of 
reference  3  and  appendix  D,  above.  "RET"  in  variable  names 
stands  for  "retarded",  thus  T.RET  is  T.RET.NXT  is 

SO  WIDTH  LPRINT  80 
110  DEFINT  I-N 
120  DIM  8(4000) 

ISO  DATA  3. 131SE10,  8.9217E-2,  1.B6988E12,  3.S0998E11,  3, 
4.18879,  339.292,  4.96E10,  3. 141S9,  .519739 
160  READ  Cl,  C2,  CQ,  CE,  CL,  C4PI3,  C4PC3,  Q,  PI,  CT 
170  HOME 

:  VTAB  2 

180  PRINT  "Enter  parameters  for  triangular  yield  function:" 


E  ~  1 


PRINT 


190  HTAB  10 

:  INPUT  “Pulse  Width  <shakes>:  T. pulse  -  ",T. PULSE 
200  T.PK  -  T. PULSE  /  2 
220  PRINT 
240  PRINT 

:  PRINT  “Enter  Burst  Parameters:" 

:  PRINT 
250  HTAB  10 

:  INPUT  “Yield  (kiloton)  ■  ",  YIELD 

:  YIELD  -  YIELD  «  4.2E-i-12  'Convert  kilotons  to  Joules 
270  YDOT.PK  ■  YIELD  /  T.PK 
280  PRINT 

:  PRINT  “Enter  Quadrature  Parameters: “ 

:  PRINT 
290  HTAB  10 

:  INPUT  "Time  Increment  (shakes)  ■  DELTA. T# 

:  DT  -  DELTA. Ttt 
300  PRINT 

:  C. DECAY  -  EXP(-C2*DT) 

310  HTAB  10 

■  INPUT  “Number  o-f  steps  per  printout  ~  N. PRINT 
320  PRINT 

:  PRINT  "Commencing  computations" 

330  SOTO  500 
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350  *Data  Lprint  Subroutina 

360  LPRINT  T,  T.RET,  R,  V,  Y. INTER 

365  T.EFF  -  CT  «  SaR(SQR(E)> 

370  LPRINT  T.EFF,E,S,U,R.OaT 
380  LPRINT  F,6,H 
390  LPRINT 
400  RETURN 

500  ’Initialization  of  variables 

510  R.DOT  -  CL 

520  R«0 

530  V-0 

550  S-Ci 

560  8(0) -8 

580  E-0 

590  1-0 

600  T-0 

610  T.RET-0 

620  IF  DT  >  T.PK 

THEN  PRINT  “Error:  pulse  shorter  than  increment 
:  BEEP  25, 100:8T0P 

EL8E  U  -  YDOT.PK  *  DT  /  2  /  T.PK 
630  F-0 
640  8-0 
650  H-0 
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h60  Y. XNTEG-0 


700  'Initial  printout 

710  LPRINT  "Burnout  Phaso  Firaball  Growth  Program:" 

711  LPRINT 

"Version  2.  Itr  —  Single  Precision,  No  Ground  Ef-fects" 
:  LPRINT  "Yield  <KT)  »  “I YIELD/4. 2E+12 

712  LPRINT 

:  LPRINT  "Pulse  Shape  Parameters:" 

713  LPRINT  "  Triangular  Pulse:" 

:  LPRINT  "  start  at  t  «  O,  peak  at  t  -"IT.PKI 
"  stop  at  t  -"IT. PULSE 

714  LPRINT  "  Peak  Yield  Rate  (J/shake)  -"lYOOT.PK 

715  LPRINT 

716  LPRINT  "Quadrature  Parameters: " 

717  LPRINT  "  Time  Increment  (shake)  -  "IDT 

71B  LPRINT  "  Number  per  printout  -  "IN. PRINT 

719  LPRINT:  LPRINT:  LPRINT 

720  LPRINT  "T  (shake) ", "T. retarded", "Radius  (m)”, 

"V  (0-^3)  ",  "Y.integ  (J)" 

730  LPRINT  "T.eff  (eV)","E-th  (J/m''3)","S  i J/m''3/sh) ", 

"U  (J/shake)", "dR/dt  (m/sh)” 

740  LPRINT  "F  (J/m''3)","G  (J-sh/m-^S)  ",  "H  ( J-sh''2/m''3) " 
7S0  LPRINT 

760  GOSUB  3S0  'data  1 print 
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1000  'Main  Iteration  Routine 

1010  R.NXT  «  R  +  R.DOT  «  DT 

1020  V.NXT  «  C4PI3  «  R.NXT  «  R.NXT  *  R.NXT 

1040  ETA  -  (V  /  V.NXT)  *  C. DECAY 

1050  E.NXT  -  ETA  *  E  +  Cl  *  DT  *  (1  ETA)  /  2 

1060  T.NXT  -  <I+1)  *  DELTA. T# 

1070  T.RET.NXT  -  T.NXT  -  R.NXT/CL 
1080  DT.RET  -  T.RET.NXT  -  T.RET 
1090  S.NXT  -  C2  *  (CE  -  E.NXT) 

IlOO  Bd-i-l)  «  S.NXT 

1110  l.RET.NXT  -  FIX (T.RET.NXT  /  DT) 

1120  8.RET.NXT  -  8 (l.RET.NXT) 

+  (T.RET.NXT  /  DT  “  l.RET.NXT) 

*  (S(I.RET.NXT+1)  -  S(l.RET.NXT)) 

1140  'Recursive  system  to  eval  integrals  o-f  S(t) 

1150  A1  -  DT  «  8  /  2 

1160  A2  «  DT  »  S.NXT  /  2 

1170  A3  -  DT.RET  *  8. RET  /  2 

1180  A4  -  DT.RET  *  8.RET.NXT  /  2 

1190  F.NXT  -F+A1+A2-A3-A4 

1200  A6  -  R.NXT  /  CL 

1210  AS  -  A6  ^  DT 

1220  6.NXT  -  8  >  DT.RET  *  (F-A3)  +  A1  *  A5  -t-  A2  «  A6 
1230  H.NXT  -  H  +  DT.RET  *  (-2*8  +  DT.RET  *  (F-2*A3/3) ) 

*  A1*A5«AS  *  A2«A6«A6 


r 

f 

1232  IF  T.RET.NXT  -  0 

THEN  YDOT.RET.NXT  -  U 

:  GOTO  1240  'Not*:  the  GOTO  is  part  of  the  THEN 

1233  IF  T.RET.NXT  >«  T. PULSE 

THEN  YDOT.RET.NXT  -  O 
ELSE  IF  T.RET.NXT  <»  T.PK 

THEN  YDOT.RET.NXT  -  YDOT.PK  *  T.RET.NXT  /  T.PK 
ELSE  YDOT.RET.NXT  »  YDOT.PK 

*  (T.PULSE-T.RET.NXT)  /  T.PK 
1240  U.NXT  -  YDOT.RET.NXT  -  C4PC3  «  CSNG(H.NXT) 

1230  R.DOT  «  CL  /  (1  CQ  «  R.NXT  *  R.NXT  /  U.NXT) 

1260  R-R.NXT 

1270  V-V.NXT 

1290  EaE.NXT 

1300  T-T.NXT 

1310  T.RET-T.RET.NXT 

1320  S-S.NXT 

1330  S.RET-S.RET.NXT 

1330  F-F.NXT 

1360  G-G.NXT 

1370  H-H.NXT 

1380  U-U.NXT 

1390  I-I+l 
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1400  *End-o-f -Phase  Criteria  checks 
1410  IF  T>T. PULSE 

THEN  Y.INTEG  -  YIELD 
ELSE  IF  T<-T.PK 

THEN  Y.INTEG  =  YDOT.PK  *T*T/2/T.PK 
ELSE  Y.INTEG  -  YIELD 

-  YDOT.PK  *  <T.PULSE-T>  *  <T.PULSE-T)  /  2  /  T.PK 
1430  IF  U<-0  OR  V*<Q+E)>»Y.INTEG 
THEN  GOSUB  350:  GOTO  2000 

1500  'Printout  Criteria  Check 
1510  IF  (I  MOD  N. PRINT)  «  O 
THEN  GOSUB  350 

1600  'Loop  to  next  iteration 
1610  GOTO  1000 

\ 

2000  'Problem  completed. 

Loop  back  to  input  routine  for  next  problem 
2010  LPRINT  CHR4<12)  'Form  Feed 
2020  BEEP  40,  250 
:  GOTO  170 

'Beep  sounds  speaker  to  notify  operator  ready  for  next  case 
2030  END  'Not  Executed,  present  for  compiler  only 

Notes  at  end  of  Appendix  A  are  applicable  to  this  program. 
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Appendix  F 

l§3i.£  Computer  Program 
■for 

Diffusion  Phase  of  the  Earl_y  Ti^me  Model 

This  version  of  the  program  assumes  a  triangular  yield 
pulse  of  isosceles  shape.  It  assumes  that  the  end  of  the 
burnout  phase  occurs  after  the  peak  of  the  yield  pulse  but 
before  the  end  of  the  pulse.  If  the  end  of  the  burnout  phase 
is  in  fact  after  the  end  of  the  yield  pulse  for  a  particular 
problem,  a  revised  version  of  the  program  (Mhich  eliminates 
the  terms  for  yield  rate  contribution  to  E^)  should  be  used. 
This  program  implements  Pomraning's  finite  differencing 
scheme  CRef  31. 

50  MIOTH  LPRINT  80 
60  DEFINT  I-N 

70  DEF  FNC0EF<T)  -  RHO. KAPPA. C 

*  (T.B/T.M)  *  (T.B/T.M)  ♦  (T.B/T.M) 

75  DEF  FNGiT)  -SO  *  (T/T.B)  *  (T/T.B)  *  (T/T.B) 

90  DATA  8.9217E-2,  13.7,  400,  3,  3.14159,  4.96E10, 

1.24E8,  296.973 

100  READ  RHO. KAPPA. C,  A,  T.B,  C,  PI,  E.B,  C.V,  60 
110  HOME 

:  VTAB  2 

120  PRINT  "Enter  parameters  for  triangular  yield  function:" 


F  -  1 


PRINT 


160  PRINT 

:  PRINT  "Enter  Burst  Parameters: " 

:  PRINT 
170  HTAB  10 

:  INPUT  “Pulse  Width  (shakes):  T. PULSE  «  ",T. PULSE 
180  T.PK  -  T. PULSE  /  2 
190  HTAB  10 

:  INPUT  "Yield  (kiloton)  »  ",  YIELD 

:  YIELD  >■  YIELD  *  A.2E+12  'Convert  kilotons  to  Joules 
193  YDOT.PK  -  YIELD  /  T.PK 
200  PRINT 

:  PRINT  "Enter  Quadrature  Parameters: " 

:  PRINT 
210  HTAB  10 

:  INPUT  "Time  Increment  (shakes)  ■  ",  DELTA. T# 

:  DT  «  DELTA. T# 

220  PRINT 
230  HTAB  10 

:  INPUT  "Number  o-f  steps  per  printout  “  N. PRINT 
240  PRINT 

:  HTAB  10 

:  INPUT  "Initial  Time:  ",T# 

:  T  -  T# 

242  PRINT 
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HTAB  10 


:  INPUT  “Initial  Counter:  “,I 
244  PRINT 

:  HTAB  10 

:  INPUT  “Maxifflum  time  to  which  to  compute:  " 
250  PRINT 

:  HTAB  10 

:  INPUT  “Yield  delivered  so  far:  ”,YSUM 
260  PRINT 

:  HTAB  10 

:  INPUT  “Initial  thermal  energy  density 
:  E.T  -  “,E.T 

265  PRINT 

:  HTAB  10 

:  INPUT  "Initial  material  energy  density 
(O  for  default):  E.M  »  ",E.M 

270  PRINT 

:  HTAB  10 

;  INPUT  “Initial  Volume:  V  -  ",V 
280  GOTO  400 

300  'Data  Lprint  Subroutine 

310  LPRINT  T,  V,  R,  R.DOT 

320  T.EFF  -  SQR  (SQR  (E.T  /  A)  > 

330  LPRINT  E.8,  E.M,  E.T,  T.M,  T.EFF 
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,T.QUIT 


mmmm 


! 


340  LPRINT 
350  RETURN 

400  'Initialization  of  variables 
410  IF  E.M  -  O 

THEN  E.M  -  4.96E+10 
420  E.S  ■  YSUM  /  V  -  E.M  -  E.T 
430  T.M  *  E.M  /  C.V 
440  R  ■  <.75*V/PI)'^<l/3) 

450  R.DOT  «  0 
460  IF  T  >  T. PULSE 

THEN  PRINT  “Error:  use  version 

’DBOIPTI'  if  pulse  is  already  over” 

:  BEEP  25,200 
:  STOP 

ELSE  IF  T<T.PK 

THEN  PRINT  "Error:  initial  time  <  t.peak 
not  supported" 

:  BEEP  25,200 
:  STOP 

465  YSUM  =  YIELD 

-  YDOT.PK  *  (T.PULSE-T>  ♦  (T.PULSE-T)  /  2  /  T.PK 
470  Q.M  -  E.M  *  V 
480  Q.S  •  E.S  *  V 
490  Q.T  -  E.T  ♦  V 
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600  'Initial  printout 

610  U»R1NT  "Diffusion  Phass  Radiativs  Fireball  Growth 
Program  based  on  Pomraning's  model:" 

620  LPRINT  "Version  1.2  —  Single  Precision 

—  Radii  assume  spherical  fireball" 

:  LPRINT  " —  Triangular  pulse  yield  integration" 
;  LPRINT  "YIELD  <KT)  «  “I YIELD/4. 2E+12 
630  LPRINT 

:  LPRINT  "Pulse  Shape  Parameters:" 

640  LPRINT  "  Triangular  Pulse:" 

6S0  LPRINT  "  start  at  t  -  O  peak  at  t  «";T.PK| 

"  stop  at  t  -“IT. PULSE 

660  LPRINT 

670  LPRINT  "Quadrature  Parameters: " 

6B0  LPRINT  "  Time  Increment  (shake)  -  "|DT 

690  LPRINT  "  Number  per  printout  -  "IN. PRINT 

700  LPRINT:  LPRINT:  LPRINT 

710  LPRINT  "T  (shake)  ","V  (m-^S) ",  "Radius  (m)", 

"dR/dt  (m/shake) " 

720  LPRINT  "E.s  ( J/m^3)  ",  "E.m  (J/m^3)  ",  "E. th  (J/m-'3)", 
"T.m  (eV)","T.eff  (eV) " 

730  LPRINT 

740  60SUB  300  'data  Iprint 
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1000  *Main  Iteration  Routine 


1010  'Compute  coe-f-f icients 
1020  IF  T.M  -  T.B 

THEN  RHO. KAPPA. C. OF. T  -  RHO. KAPPA. C 
:  G  -  60 

ELSE  RHO. KAPPA. C. OF. T  -  FNCOEF(T.M) 

:  6  -  FN6<T.n> 

1030  AO  -  RHO. KAPPA. C. OF. T 
1040  IF  T+DT  >  T. PULSE 

THEN  YSUM.NXT  -  YIELD 
ELSE  YSUM.NXT  >■  YIELD 

-  YDOT.PK  *  (T.PULSE-T-DT) 

♦  (T.PULSE-T-DT>  /  2  /  T.PK 

1050  DY  -  YSUM.NXT~YSUM 
1060  YSUM  «  YSUM.NXT 

1070  EPS  -  <G*a.T*DT/V'^<2/3)  -  DY>  /  YSUM 
lOSO  IF  EPS  <  .001 

THEN  V.NXT  -  V  *  <1  +  3  ♦  EPS) 

ELSE  V.NXT  -  V  *  <1  +  EPS)  ♦  (1  +  EPS)  *  (1  ♦  EPS) 
1090  A1  «  1  <«•  AO  *  DT 
1100  A2  «  A  «  T.M  ♦  T.M  *  T.M  /  C.V 

1200  'Compute  new  energy  partition  among 

streaming,  material,  thermal  radiation 
1210  Q.S.NXT  -  (Q.S  DY)  /  A1 
1220  A3-14-A0«DT«  <1-»'A2) 

F  -  6 


1230  Q.M.NXT  •  (Q.M  *  A1  Q.S.NXT  «  AO  «  OT  «  A1 

+  Q.T  «  AO  *  DT)  /  A3 

1240  Q.T.NXT  -  (Q.T  ^  AO  *  A2  *  Q.M.NXT  «  DT>  /  A1 

1300  *S«t  variables  to  next  values 

1310  Q.M  -  Q.M.NXT 

1320  Q.T  -  Q.T.NXT 

1330  Q.8  «  Q.S.NXT 

1340  E.M  «  Q.M  /  V 

1350  E.T  -  Q.T  /  V 

1360  E.8  -  Q.S  /  V 

1370  T.M  -  E.M  /  C.V 

1390  V  -  V.NXT 

1400  R.NXT  -  (.75  *  V  /  PI)  (1/3) 

1410  R.DOT  «  (R.NXT  -  R)  /  DT 

1420  R  -  R.NXT 

1430  T#  -  T#  +  DELTA. T* 

1440  T  -  T# 

1450  1  -  1  1 

1500  'End-of-Phase  Criteria  checks 
1510  IF  T.M  <-  25  OR  T>-T.QUIT 
THEN  60SUB  300 

:  60T0  1600  *Note:  GOTO  is  part  of  THEN 
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1S20  ’Printout  Criteria  Check  ' 

1S30  IF  (I  MOD  N. PRINT)  «  O 
T)£N  I  »  O 

:  60SUB  300 

1S40  'Loop  to  next  iteration 
1550  GOTO  1000 

1600  ’End  of  Radiative  only  -fireball  growth.  * 

Hydroseparation  assumed  to  occur  at  temp  of  25  eV 
or  T.quit  was  reached  first 
1610  LPRINT  CHR«(12)  ’Form  Feed 
1615  BEEP  40,250 

:  GOTO  110  ’Loop  back  to  input  routine  for  next  case 
1620  END  ’Not  excK:uted,  present  for  compiler  only 

Notes  at  end  of  Appendix  A  are  applicable  to  this  program. 
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